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FOREWORD 


Backg rou,.d - The RF technique is a gauging method which 
samples all parts of the inside of a tank. The response is characteristic 
of the total mrss of fluic! within a tank; this is one of the advantages cf 
the RF system over the standard capacitance system which is essentially 
a local measurement of fluid density. Another advantage of the RF sys- 
tem is the simplicity of the hardware involved; a small grounded antenna 
about the size of a paper clip is sufficient to communicate with the inside 
of the tank. 

Preliminary theoretical and experimental results on the RF gauging 
idea were obtained by NBS in connection with a NASA sponsored contract 
on slush hydrogen gauging. These results may be found in NBS Report 
9793 dated June 1, 1971, on "Instrumentation for Hydrogen Slush Storage 
Containers. " 

Purpose - The purpose of this report is to summarize work done 
under purchase order T-1738B from the NASA Johnson Space Center 
Houston, Texas, to the National Bureau of Standards, Cryogenics Division, 
Boulder, Colorado. Items covered includes 

1) Phase I - Preliminary studies of the radio frequency (RF) 
mass quantity gauging system for two phase and supercritical 
fluids; and construction of experimental system for detailed 
feasibility studies. 

2) Phase II - Experimental evaluation of the system for super- 
critical nitrogen and hydrogen. (Oxygen is also included in 
Phase II; the results will be reported separately upon the 
completion of the oxygen testing. ) 

PEECffiING PAGR bunk not films, 
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Objective - The primary objective of this work is to design and 
develop a breadboard system to verify that the radio frequency resonant 
cavity mode analysis technique is conceptually sound for the fluid mass 
quantity gauging of the Space Shuttle Orbiter PRSD (Power Reactant 
Storage and Distribution) subsystem tankage, i. e. , supercritically stored 
hydrogen and oxygen, in all gravity fields. The secondary program 
objective is to analytically determine the applicability of the concept to 
the quantity gauging of Shuttle Orbiter propulsion systems tankage. 

(i. e. sub-critical fluids). 

End Product - The end product of this contractual effort is to be 
a breadboard RF mass quantity gauging system capable of gauging super- 
critically stored hydrogen and oxygen to an accuracy of one percent of 
total tank quantity in any gravity environment. This should also include 
the gauging of subcritical hydrogen and oxygen representative of fill and 
drain operations on the PRSD tanks. 

Acknowledgements - The authors wish to acknowledge D. B. Mann, 
Chief of the Cryogenic Metrology Section, Peter V. Tryon and John F. 
LaBecque of the NBS Statistical Engineering Laboratory, and Dan S. Trent, 
the NASA, JSC technical monitor, for assistance and encouragement in 
planning this program. The authors also appreciate the assistance of 
Larry Lowe, Bill German, Bob Gray, Tim Hipsher, Ed Rogers, Brian 
Nakagawa, and Jeannette Garing. 
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ABSTRACT 


This is a summary report of work done to date on NASA (Johnson 
Space Center) purchase order T-1738B concerning Radio Frequency (RF) 
Mass Quantity Gauging. Experimental apparatus has been designed and 
tested which measures the resonant frequencies of a tank in the "time 
domain. " These frequencies correspond to the total mass of fluid with- 
in the tank. Experimental results are discussed for nitrogen and hydro- 
gen in normal gravity both in the supercritical state and also in the two 
phase (liquid-gas) region. Theoretical discussions for more general 
cases are given. 

Key Words: Gauging; hydrogen; nitrogen; radio frequency; total mass. 



INTRODUCTION 


When a small antenna is placed in a closed metal cavity, the 

electromagnetic field pattern which the antenna generates inside the 

cavity depends on the excitation frequency and the shape of the cavity. 

At certain frequencies, f (called resonant frequencies), the field 

n 

patterns are standing waves (called resonant modes). These modes are 
easily detected at the antenna since the impedance match of the antenna 
to the cavity is more efficient at the resonant frequencies. 

The presence of a dielectric fluid within the cavity will change 
the resonant frequencies. The resonant frequencies will decrease with 
an increasing amount of fluid because the presence of the fluid slows 
down the propagation of the electromagnetic wave. This presents the 
possibility of gauging the amount of fluid by measuring the resonant 
frequencies. 

If the density of fluid is uniform throughout the cavity. 


f 

n 



( 1 ) 


where f is the resonant frequency of the empty cavity for the nth 
mode and * is the dielectric constant of the fluid. For many nonpolar 
fluids of interest (including hydrogen and oxygen), the dielectric con- 
stant depends only on the density of the fluid. In this case there is a 

unique relationship between each resonant frequency and the total mass 
within the cavity; and only one mode is necessary to determine the total 


mass. 



If the density of the fluid is not uniform throughout the cavity 
(either because of a two phase liquid-gas interface or a single phase 
fluid with temperature gradients) the resonant frequencies depend on 
the amount of fluid mass in the cavity and also somewhat on where the 
dense portions of the fluid are located within the cavity (fluid geometry 
effects). Since the geometry of each standing wave is different (the 
modes are linearly independent functions of the space variables), it is 
possible to partially compensate for the fluid geometry effects by com- 
paring two or more modes. This process of comparison is called Mass 
Quantity Gauging by RF Mode Analysis. 

The purpose of this present work is to develop an experimental 
system which will 

(1) Provide a breadboard total mass gauge for measuring 
uniform density fluids; developing the accuracy and data reduction of 
the time domain technique for measuring a single resonant frequency, 
and 

(2) Determine the feasibility of RF Mode Analysis for nonuniform 
fluids; to measure the total mass with sufficient accuracy for dynamic 
Normal-g and Zero-g fluid geometries. 

This report will emphasize item (1) above under the heading 
Uniform Density Fluids and also give some preliminary discussion and 
results of item (2) above under the heading Nonuniform Density Fluids. 



EXPERIMENTAL SYSTEM 


The experimental system consists of (1) an electronic signal 
conditioner and data acquisition system which measures the resonant 
frequencies in the time domain and (2) an 18 inch diameter spherical 
cryogenic storage tank (experimental vessel) designed for testing nitrogen 
and hydrogen, in both the two-phase and supercritical states; it is ex- 
pected that after further testing that the tank will also be suitable for 
oxygen. The total fluid mass is determined by weighing with a calibrated 
load cell; there are a number of thermocouples and resistance thermom- 
eters attached to the sphere to measure temperature gradients in the 
fluid. There are several antenna locations for measuring the effect of 
antenna orientation with respect to a non-homogeneous fluid. The details 
of the cryogenic system are outlined in Appendix A. 

The resonant frequencies may be detected by sweeping the antenna 
with an RF sweep generator; where the generator frequency ranges be- 
tween f and fg and the resonant frequency (or frequencies) of interest 

lies between f and f . When the generator frequency coincides with the 
A B 

resonant frequency there is a decrease in the signal reflected from the 
antenna which shows up as a spike in the detector output; this output may 
be displayed on an oscilloscope (see figure 1). 

If the sweep generator frequency output is linear (or at least 
repeatable) in time then the resonant frequencies may be measured by 
measuring the time interval between the output of a reference cavity 
(tuned to a frequency f^) and the output of the experimental vessel (see 
figure 1). For example, if the sweep rate r is linear, the resonant 
frequency of the fundamental mode, f is given by 


V 


'*1 - * 0 ' 


( 2 ) 


3 






Since f and r are iixed, fj is determined by the time interval ( t ^ t^); 

this then is a measurement of the resonant frequency in the "time dorr'. i 1 . 

The time intervals can be measured using a digital clock and a 
counter. The clock may be triggered in a start and stop process by 
signals coming from the reference cavity and experimental vessel, re- 
spectively. This is shown schematically in Figure 2. The details of the 
signal conditioner and data acquisition system are contained in Appendix B. 

If the detectors are to start and stop the clock in a precise manner 
then the pulse must be very sharp and narrow so that the signal in real 
time comes precisely at the time the RF generator is at the resonant, 
frequency. The narrowness of the pulse is related to the Q of the cavity 
which is defined by 

f 

Q = TT • (3) 

n 6f 

n 

where 6f^ is the width of the spike at the hall power points. For 
example, if is 10, 000 then there will be about a 0. 01 percent un- 
certainty in the measurement of f . 

n 

The Q for several of the resonant modes have been measured 
n 

in detail for a 19 inch diameter copper sphere, an 18 inch diameter 
stainless steel sphere and a 5 foot diameter stainless steel sphere. The 
measured Q values range between 6,200 and 91,000. The detailed 
results are contained in Appendix C. 

The conclusion of Appendix C is that the Q's are high enough to 
accurately measure the resonant frequencies by the time domain technique 
even for the large vessel; and it is reasonable to use small vessels for 
scaling experiments on ultimate large tank configurations. 
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Figure 2. Conversion of resonant frequencies to the "time domain. " 






THE RF ANTENNA 


One of the attractive features of the RF technique is the 
simplicity of the internal tank hardware, which is simply a small loop 
of wire. Figure 3 shows two of the antenna configurations which have 
been used ir. the experimental vessel; these are connected to high pres- 
sure coaxial feedthroughs. 

The straight wire antenna is the "TM probe". It generates 
only the TM modes. The straight wire is simply an extension of the 
center conductor of the coaxial feedthrough. 

In the loop antenna, the center conductor is bent into a U- shape 
about 3/4 inch by 3/4 inch and the end is grounded to the outer con- 
ductor. This antenna will generate both TE and TM modes. 

The coupling of the antenna to the cavity (ard hence the amplitude 
of the response) is changed only slightly by changes in the size and shape 
of the antenna. There appears to be a wide variety of antenna size and 
shapes which are acceptable for this gauging technique. 
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For a non-polar dielectric fluid of density, c. (a constant through- 
out the cavity), the resonant frequency of the nth mode f is given by^ 


f - ° n 
n /T(o) 


) P 
i i 


u it 


where f is the empty cavity frequency. Semi-emperically e(o) can be 


on 


given implicity by^’ 


= A 0 + 

e( o) + 2 


Bo 2 + Cp 3 


( 5 ) 


where A, B, and C are constants determined experimentally for each 
fluid. The f are determined experimentally and serve to calibrate 
the system, f depends on the size and shape of the cavity; if the cavity 
changes size because of thermal contractions or pressure expansions f 
must be adjusted accordingly; f may also change if objects are placed 
in the cavity. To a good approximation, for most situations of interest B 
and C may be neglected in equation (o); if V is the total volume of the 
cavity, the mass M in this case may be given by 


2 ^2 

. . V / e(p)-l \ V_ / on /n \ 
M = A^(o) + 2/-A +2f 2 ) 

on n 


( 6 ) 


corrections to this formula for non zero B and C may be applied if 
necessary. It is seen that possible inaccuracies in the total mass come 
from four sources: 

1. The uncertainties in A, B, and C. This is not a serious 
problem if the properties data taken by careful capacitance measurements 




1, 2, 3, 4 See references on page 32. 
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are complete over the ranges of interest. This must be determined 
for each application involving a specific fluid. 

2. The volume of the container. This is usually inferred by 
weighing a fluid with a known density. 

3. The accuracy of f which can be measured accurately for an 

on 

empty cavity; however, if the cavity changes size or shape as a function 
of fluid density the empty cavity value is no longer valid in equation (6) 
and corrections must be made to equation 6 to account for this fact. 

4. The uncertainties in f^. For normal RF frequency ranges and 
high Q cavities this measurement can be made very accurately. The inac- 
curacies in converting this frequency into a useful digital or analog signal 
are presently about 0. 2 percent full scale; this can be improved if necessary. 

It should be noted that the factors 1, 2, and 3 above may be by- 
passed by direct calibration of f vs M using a gravimetric weigh 
system. However, this is not practical to do for every system. The 
purpose of this work, then, is to examine the validity of equation (6) (and 
possibly corrections to equation (6)) for general use in a system where 
the four accuracy factors may be adequately evaluated. This will be 
described for the cases of nitrogen and hydrogen. It is anticipated that 
oxygen will also be evaluated in the near future. The general approach 
will be to directly calibrate f vs M using a gravimetric weigh system 
for measuring M and a calibrated reference cavity for measuring f 
and comparing these results with equation (6) using the "time domain" 
method of measuring f^. The amount of any possible density variations 
will be inferred from an array of thermocouples attached to the cavity. 

Data Reduction and Readout 

Although equation (6) is fairly complicated in form it is surpris* 
ingly linear over the density range of interest. A least squares linear fit 



to equation (6) for the density range from zero to the normal boiling 
point shows a maximum deviation of 2, 45 percent of the total range for 
oxygen, 2. 2 percent for nitrogen and 1. 3 percent for hydrogen. Thus 
the frequency vs mass relation can be expressed to a good approxi- 
mation by 


f = f 
n on 


(ot + Bo) 


(7) 


where ev and 0 are constants depending only on the fluid. This simple 
linear readout may be sufficient for many applications. Quadratic fits 
of the form 

f = f (<v + 00 + VO 2 ) (8) 

n on 


may be found which give less than 1 percent error for all fluids (see 
Appendix E). 

Experimental Data for Nitrogen and Hydrogen 

Preliminary experiments were started with liquid nitrogen to 
check out the experimental apparatus and the theory of uniform density 
fluids. The pressure was raised to 500-700 psi and the experimental 
vessel was agitated to achieve nearly uniform density. The frequencies 
were measured directly using the reference cavity; the line of the reference 
cavity was shifted to coincide with the resonant line on the face of the oscil- 
loscope. Figures 4 and 5 show the data taken in this manner for the TM 0 n 
and TMo^ modes respectively. It is seen that the data are nearly linear 
(slightly concave downward) and corresponds very well to the theoretical 
form given by equation (6). We believe that the occasional point which 
deviates from the theoretical curve is due to non-homogeneous density 
"stratification. " We will give some preliminary results on stratification 
for the case of hydrogen. 


11 



TOTAL MASS, lbs 



Figure 4 # Mass gauging of uniform density LN 2 using the TM 0l i mode 








Figure 5. Mass gauging of uniform density LN ? using the TM 0?1 mode, 
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The hydrogen data were taken with the data acquisition system 
described in Appendix B. The resonant frequencies for the first five 
fundamental modes were converted to milliseconds in the time domain 
(as indicated in equation 2) to the nearest 0. 01 msec. These times 
were digitized by an internal clock and stored on magnetic tape. The 
tape was used to punch data cards which contain a run code, the time 
conversion of the resonant frequencies, and the scatter in the time domain 
data. The scatter in the time domain data (except for an occasional noise 
spike) was ± 0. 01 msec corresponding to the first significant digit. The 
data cards also contain the pressure, temperature and mass data cor- 
responding tc the run identification number. Computer plots of the data 
for the TM 01 i and TM^ modes are shown in figures 6 and 7 respectively. 
A least squares fit of these data to equation (6) gives a 3a (99. 9% confi- 
dence level) for the TM on mode of 1. 2 percent, the maximum deviation 
of the data points taken is 1, 10 percent. Further details of the data 
analysis, and also the accuracy statements are contained in Appendices 
D and E. 
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Figure 6. Mass gauging of uniform density LH # using the TM 
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Figure 7. Mesa gauging of uniform density LH, using the TM ot i mode. 
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EFFECTS OF STRATIFICATION 


In the data described above, every effort was applied to achieve 
uniform density; however, a few preliminary data have been obtained 
concerning the effects of stratification on a single resonant frequency. 
The following cases (in which the top and bottom temperatures were not 
equal) gives a rough idea of the effect of stratification on the TM on mode 
for supercritical hydrogen at 350 psi: 


Case I 


Top temperature 

= 

47. 5 K 

Bottom temperature 

= 

27. 5 K 

Top density 

A*/ 

1 lb /ft 3 

Bottom density 


4 lb/ft 3 

Average density 


3. 02 lb /ft 

Measured weight (load cell) 

= 

5. 35 lbs 

Inferred weight (from fig. 6) 

= 

5. 55 lbs 

% Full scale error 


4% 


Case II - Top temperature = 

Bottom temperature = 

Top density — 

Bottom density ~ 

Average density ~ 

Measured weight (load cell) = 
Inferred weight (from fig. 6) = 


37. 5 K 
33 K 

2. 25 lbs /ft 3 

3. 45 lbs /ft 3 

2. 00 lbs /ft 3 

3. 55 lbs 
3. 45 lbs 


% Full scale error 1. 5% 


Case III - Top temperature 

Bottom temperature 
Top density 
Bottom density 
Average density 


= 31. 5 K 
= 53 K 
~ 3. 8 lbs /ft 3 
~ 0. 8 lbs /ft 3 
- 1. 07 lbs /ft 3 
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Measured weight = 1.9 lbs 

Inferred weight (from fig. 6) = 1. 9 lbs 

% Full scale error ~ 0. 0% 

In Case I, the stratification is cold fluid on the bottom, warm on 
top with a good share of the fluid in between at the colder temperature; 
definitely not a linear thermal gradient. 

In Case II, there is a cold fluid on top, a colder fluid on the bottom 
and warmer fluid in the middle as inferred from the average density. 

In Case III, there is a cold fluid on top, a warm fluid on the 
bottom with most of the fluid in the cavity at the warmer temperature. 

From these cases, it is seen that certain small amounts of strati- 
fication may be tolerable; this is because the antenna senses the entire 
cavity and tends to have an integrative effect over all the mass within 
the oavity. For larger amounts of stratification, one antenna and one 
mode may not be sufficient to achieve the desired accuracy, and more 
information from other antennas or modes may be necessary to properly 
use the RF technique as a gauge. This situation is discussed in the next 
section on non-uniform density fluids. 


I 
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NON-UNIFORM DENSITY FLUIDS - RF MODE ANALYSIS 


Non-uniform (inhomogeneous) density (or dielectric constant) may 
occur either in a single phase fluid under temperature gradients or in a 
two phase fluid with a liquid-gas interface. The spacial gradients in the 
dielectric constant have a diffractive effect on the propagation of the 
electro-magnetic wave and this changes the shape of the standing wave 
patterns of each resonant mode. These change:' will be different for each 
mode because of the dissimilar field patterns of the modes. The resonant 
frequency of the nth mode will be given by 

f = f + Af 
n un n 

Where f is the resonant frequency expected for the uniform density case 
assuming that the total r.. a ss in the cavity is spread uniformly over the 
cavity; Af^ is the change in resonant frequency due to the non-uniform 
geometry of the fluid and will be different for each mode even to the extent 
of being positive or negative depending on the mode. 

It is natural to ask: what are the extreme limits of Af as the fluid 

n 

ranges over all possible configurations? This is a difficult and possibly 

impractical question to answer. The reason is that for a given mass of 

fluid, there are two fluid geometries which will give the maximum and 

minimum values for Af l hut these geometries appear to be very compli- 

n r 

cated and it is unlikely that they will occur in practice It is more feas- 
ible to talk about the practical limits of Af which are determined by 

n 

experimentation and calculation of A f for likely fluid geometries. 

n 

Theoretical Approach 

Theoretical work which was initiated in early stages of Phase I 
was directed along two lines. The first was to develop approximation 
techniques for calculating how fluid geometry affects the resonant 
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frequencies; the second was to investigate the mode geometries to get a 
qualitative picture of how fluid location may affect the resonant fre- 
quencies. 

Earlier work had calculated in closed form the resonant fre- 
quencies expected for a two phase fluid with a concentric spherical phase 
boundary ("zero-g" geometry). This work is contained in Appendix F 
for reference purposes. Further work on two phase geometry effects 
must be handled by approximation or numerical techniques. Appendix G 
is a survey of relevant approximation techniques and how they may be 
applied to the cavity problem. Several of the examples are worked out 
for the case of the normal two phase fill geometry. One of the conclusions 
of this work is that the resonant frequencies are most affected when the 
dense portion of the fluid moves in and out of the high field region. These 
high field regions are different for different modes. The field profiles 
for a few of the lowest order modes are plotted in figures 8 through 13. 

It is seen that the modes partition the cavity into several diotinct high 
field regions; from this, it is expected qualitatively that an average of 
the resonant frequencies of the lowest order modes may give a mass 
value that is relatively independent of the location of the dense fluid. 

The numerical derivation of these graphs as well as numerical solutions 
of some of the approximation techniques are contained in Appendix H. 

Experimental Approach 

The experimental system described earlier is designed to record 
on magnetic tape the resonant frequencies of the first five modes for 
three different antenna locations; the data for each antenna is collected 
every 0. 1 sec. thus making it possible to study dynamic effects where 
the fluid is in motion. It is anticipated that this system will give useful 
information in "zero-g" simulation experiments. Some preliminary data 
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has been taken for normal gravity LN two phase fill. Since the liquid 
surface breaks the spherical symmetry, it is found that some of the 
resonant lines split into two or more closely spaced lines; some of the 
modes do not split. For example, the TM on mode stays as a single line 
and the frequency shifts as a function of total mass (see figure 14) giving 
some idea of the magnitude of Af^. Figure 15 shows the response of the 
TMoa^ mode which splits into three lines during normal fill. In this case, 
a straight forward average of the three modes is very close to uniform 
density curve for this mode as shown in figure 16. 

These two cases indicate two alternate methods of gauging the 
situation of normal fill together with uniform density: 

(1) In the case of TM 011 mode, the readout can be designed for 
uniform density and then a correction factor can be applied for the nor- 
mal fill condition. 

(2) In the case of the TM 0? i mode, the readout can be designed 
for uniform density, with electronic averaging of the three split resonant 
lines in the TM 0al time frame. 

It is reasonable to expect that both of these techniques could be develop- 
ed to give readout accuracies on the order of 1 percent; however, the 
averaging technique may be more useful if it can be generalized to tilt 
geometries and "low gravity" geometries. 
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Figure 15. Normal fill data for liquid nitrogen TM, 
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CONCLUSIONS 


We have demonstrated using hydrogen and nitrogen that the 
RF Mode Analysis Technique is conceptually sound for uniform density 
fluids and that there are encouraging results for non-uniform density 
fluids. The results for uniform density fluids should apply to any size or 
shape of tank as long as the factors listed on pages 8 and 9 can be ade- 
quately evaluated. Further work should be done on non-uniform density 
fluids for the spherical tank and also on tank shapes which are not 
spherical. 
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APPENDIX A 


CRYOGENIC SYSTEMS 

The basic system includes a pressure vessel, a heater, a pi; ng 
system and an insulation system. These components were common ;o 
all the test configurations. The normal gravity calibrations utilized a 
load cell to measure the mass accurately. The design of th flight pallet 
for the zero gravity tests included an additional storage dewar. The 
major components of these systems are described. 

The Pressure Vessel 

The vessel is a sphere constructed of 304 stainless steel. It 
was formed by welding together two spun hemispherical heads. The inside 
diameter is 17 - 1 3/ 1 6 inches. It was designed and constructed per 
the A.S.M.E. Code for Unfired Pressure Vessels. The minimum wall 
thickness is 0.21 inches and the maximum allowable working pressure 
is 750 psig as defined in this code. The pressure vessel has been hydro- 
statically tested to 1-1/2 times its working pressure at room temperature. 
Access ports for instrumentation and piping are provided. All ports, 
except the port on the bottom, are 3/4" female pipe threads. The bottom 
port is 1-1/2" female pipe thread to accommodate the heater. All ports 
are reinforced with welding spuds. 

Tfrfi Mittgrg 

Two heater configurations have been used. The first was designed 
for use with liquid nitrogen and the second for liquid hydrogen. Both were 
designed to pressurize the fluids to their critical pressures in forty 
minutes. Both were constructed with temperature independent resistance 
wire. Pressurization lines as a function of heater power and current are 
presented in figure Al. 
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Figure Al. Heater power and current vs, pressurizati 




The nitrogen heater syste , consisted of two separate one kilo- 
watt heaters that could be separately controlled. Each heater was 
wound with ten feet of 22 gage wire. Each had a room temperature 
resistance of 12 ohms. The wire was wound on two si pa rate toils. 

The coils were then inserted in the bottom access port with the wire 
in direct contact with the nitrogen. 

The hydrogen heater was constructed differently because of 
safety considerations with the flammable fluid. This heater consisted 
of one 300 watt coil of wire. It was wound with four feet of 28 gage 
wire with a room temperature resistance of 20 ohms. The wire was 
covered with ceramic beads and inserted into a coil of 1/4” O.D. 
copper tubing. The tubing was then soldered on the outside of a brass 
spindle which was screwed into the bottom port. The heater element 
does not contact the hydrogen but heats the pressure vessel from the 
outside. The copper tubing is purged wi* l helium gas and is electri- 
cally grounded. 

The Normal Gravity Calibration Tests 

For these tests, the pressure vessel was hung in a cubical frame 
constructed of aluminum angle. The piping system and the insulation 
system was also supported from the frame. The frame, and the systems 
it supported, was suspended from the load cell. All interface piping 
and instrument lines were flexible to minimize their effect on the load 
cell. Figure A-2 is a photograph of this system without the insulation 

The piping was designed for use with liquid hydrogen and was 
used with liquid nitrogen as well. The system is equipped with a 750 psig 




relief valve. The system has been pressure tested to 700 psig with 
liquid nitrogen as the test fluid and it has been tested to ^ 00 psig with 
liquid hydrogen . 

The insulation system was made up of six inches of fiberglass 
batting. To prevent cryopumping the fiberglass was covered and 
sealed with aluminized mylar. The insulation was then purged with 
helium gas and a slight overpressure was maintained during the tests. 

The load cell used was a strain gage type rated at 300 pounds 
maximum load; it has a load sensitivity of 0. 1 mv/lb with a 10 volt 
excitation voltage. Using a high precision digital voltmeter on the 
output, the resolution is approximately 0. 01 lbs. 

The pressure vessel was instrumented with several temperature 
measuring transducers to determine the temperature and density gradients 
during the supercritical tests. Nine copper- constantan thermocouples 
were attached to the outside wall of the pressure vessel at six inch 
vertical increments. Platinum resistance thermometers were inserted 
into the fluid at the top and bottom to determine actual fluid temperature. 

The Zero Gravity Tests 

The design of this system utilizes the same pressure vessel, heater 
and support frame as the normal gravity system. The piping and insulation 
systems have been changed. 

The piping system flow diagram is shown in figure A3. All the 
valves are manually operated. The design utilizes a common overboard 
vent and liquid dump, which will he compatible with the aircraft. The 
system can be refilled and pressurized in flight if necessary. 




Figure A3. Cryogenic flow system for zero-g simulation tests. 
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The complete system will be packaged in an aluminum paneled 
container six feet long by three feet high and three feet wide. The 
total weight will be approximately 700 pounds. 

The system is designed to withstand an acceleration load of 
16 g’s in any direction. This will be accomplished by using a polyure- 
thane foam insulation system. After the entire system has been assem- 
bled in the paneled container, the foam will be poured, filling all re- 
maining space. The foam, then, will not only insulate but will support 
the system. The outside paneling shall be 3/8 inch aluminum plate. 

With this insulation system, the total heat leak will be approxi- 
mately 15 BTU per hour. The rate of pressure rise, as a result of 
this heat leak, will be 1. 5 psi per hour. 
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APPENDIX B 


MASS TANK GAGING SIGNAL CONDITIONER AND DATA 
ACQUISITION SYSTEM 

Purpose 

To condition the resonant frequencies of the mass tank, to 
measure these frequencies as a function of time and store the data in a 
magnetic tape unit. 


Operation 


A radio frequency generator is swept through a designated fre- 
quency spectrum of a starting frequency f , to an ending frequency f . 

A B 


This generator is connected to a cryogenic mass tank. The mass 
tank has resonant frequencies (modes) which are enhanced each time the 
tank is energized with the generator sweep of selected frequencies. 


The design objective of the data acquisition system and signal 
conditioner is to measure these enhanced resonant modes with respect 
to a known frequency and record them. 


In the system this measurement is accomplished by using a time 
measurement technique. The generator sweep is used as a time base for 
all of these measurements. 


A reference cavity tuned to a frequency, f^, which is lower than 
the lowest cryogenic cavity mode, is the reference starting time for the 
time measurement of the mass tank system. 

When the generator energizes the reference cavity tuned to f Q , a 
pulse is generated by the cavity. This pulse is used to start a counting 
sequence using a 100 kHz clock as the basic counting device. 

The first resonant mode, f^, from the mass tank strobes the 
counter into a buffer register which loads the data to a shift register. 
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The shift register then sends the data to a magnetic tape unit. All these 
data transfers occur before the next resonant mode of the mass tank is 
generated. 

The counter continues counting with its data being strobed into 
the shift register every time a resonant mode is generated. 

Upon completion of the generator sweep, the counter is reset, 
another antenna in a different location on the mass tank is multiplexed, 
and the generator initiates a new sweep. Beginning of the sweep initiates 
a new set of data. 

Using this technique, a measurement of mass in the cryogenic tank 
is obtained as a function of time with respect to the reference cavity tuned 



Data Acquisition Schematic I 


Functional Description 

This schematic shows the timing sequences necessary to start 
the counting sequence and data transfer functions. All timing events are 
referred to the initiating pulse, t^, generated by the reference cavity. 

The pulse t Q , opens the gate to the counter and enables the real 
time clock to send 100 KHz pulses to the six stage counter. 

The next timing event occurs when the cryogenic mode pulse t 
is generated. This pulse generates a strobe command to the storage 
buffer, a parallel load command to the shift register, a serial-parallel 
mode control pulse to the shift register, and a shift command gate enable 
pulse to the serial shift, clock. 

The mode counter generates a clock disable signal at the end of 
the last mode to be measured and resets itself. 
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Data Acquisition Schematic II 


Functional Description 

This schematic shows the timing sequence necessary for control 
of the mode counter, antenna multiplexer, sequence counter and magnetic 
tape recording unit. All timing events are referred to the mass tank 
pulse generated by the modes f . 

A recording sequence is initiated by depressing the manual reset 
and manual start push buttons. The manual reset button resets all counters 
and controls flip flops to initial status for recording and controlling. The 
start push button starts the magnetic tape unit. 

The t pulses are counted by the mode counter until a preselected 
output is generated in the counter. This output multiplexes the antenna 
or sends a counting pulse to the sequence counter if the antenna multi- 
plexer is disabled. 

The antenna multiplexer generates an output for every complete 
cycle of antenna multiplexing. This output sends a counting pulse to the 
sequence counter when multiplexing is enabled. 

The sequence counter has a preselected number of counts it can 
accept before it overflows and stops the magnetic tape recorder. 
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COUNTER-SHIFT REGISTER CARD A 


F rom 

To 


Function 

i ondition 

Ai 

AM.-C3 


load shut register 

i nput 

A 5 

Dl 5-C4 


8l s serial shift register 

input 

A 2 

132 


♦ s \ olt di s t r i but ion 


A4 

Tons thumbwhee 1 

switch 

date- 1 

nput 

At 

Tons thumbwheel 

swiii h 

date- 2 

i nput 

Ab 

lens thumbwheel 

swiU h 

date- 3 

t nput 

A l 0 

Fens thumbwheel 

switch 

dati - H 

r ■ ' 

A \ 2 

Ones thumbwheel 

switch 

da f e- l 

np*- • 

A 1 4 

Ones thumbwheel 

switc h 

date- 2 

input 

Alt 

Ones thumbwheel 

switch 

date- 4 

1 T.pUt 

Alb 

Ones thumbwheel 

switch 

date - 8 

nput 

A 20 

Tens thumbwheel 

switch 

run no. 1 

input 

A22 

Fens thumbwheel 

switch 

run no. 2 

input 

A24 

Tens thumbwheel 

switch 

run no. 4 

i nput 

A26 

Tens thumbwheel 

switch 

run no. 8 

i nput 

A28 

Ones thumbwheel 

switch 

run no. 1 

input 

A 30 

Ones thumbwheel 

switch 

run no. 2 

in; ut 

A 32 

Ones thumbwheel 

switch 

run no. 4 

nput 

A 34 

Ones thumbwheel 

switch 

run no. 8 

input 

A 3t 

C26 


antenna identification 1 

input 

A38 

C28 


antenna identification 2 

input 

A40 

A42-AO0 


antenna identification 4 

ground input 

A42 

A40-A60 


antenna identification 8 

ground input 

A44 

D44 


shift register output 1 


A4b 

D46 


shift register outpu; 2 


A48 

D40 


shift register output 4 


A50 

T42 


shift register output 8 


A52 

C 52 


2^iS shift register mode control 

* r.put 

AS4 

D48 


real time counter reset 

i nput 

A5b 

A3-C3-G44 


load buffer register 

input 

A58 

D58-G56 


gated oscillator 100 K Hz clock 

input 

AtO 

A40-A42- B60 

common distribution 
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DATA CONVERTER-CONTROLLER CARD M B" AND CARD "D" 


From 

To 

Function 

Condition 

B2 

A2-C2 

+ 5 volt distribution 


B30 

F47 

magnetic tape 1 

output 

B32 

F4 C 

magnetic tape 2 

output 

B34 

F43 

magnetic tape 4 

output 

B36 

F41 

magnetic tape 8 

output 

B38 

F39 

magnetic tape A 

output 

B40 

F37 

magnetic tape B 

output 

B42 

D56 

100 K Hz clock free running 

output 

B44 

D32 

shift register data 1 

input 

B46 

D34 

shift register data 2 

input 

B48 

D36 

shift register data 4 

input 

B50 

D38 

shift register data 8 

input 

B54 

356-B60 

code conversion enable 

input 

B60 

A60-C60 

common distribution 


D28 

C40 

dummy data enable 

input 

D30 

C42 

real data enable 

input 

D32 

B44 

shift register data 1 

output 

D34 

B46 

shift register data 2 

output 

D36 

B48 

shift register data 4 

output 

D38 

B50 

shift register data 8 

output 

D40 

A48 

shift register data 1 

input 

D42 

A50 

shift register data 2 

input 

D44 

A44 

shift register data 4 

input 

D46 

A46 

shift register data 8 

input 

D52 

G24 

reference cavity 

input 

D56 

B42 

100 K Hz clock free running 

input 

D58 

A58-G56 

100 K Hz clock gated 

output 
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B1 

9311 

one of sixteen decoder 

B2 

7430 

8- input positive NAND gate 

B3 

7430 

8-input positive NAND gate 

B4 

7430 

8-input positive NAND gate 

B5 

7420 

dual 4-input positive NAND gate 

B6 

7400 

quad Z-input positive NAND gate 

B7 

2 capacitor cambion 


B8 

832 

dual buffer 


D17 

2 K resistors 



D21 

7404 

HEX inverter 


D22 

7403 

quad 2-input positive 

NAND gate 

D23 

7403 

quad 2- input positive 

NAND gate 

D29 

848 

clocked flip-flop 


D30 

7400 

quad 2- input positive 

NAND gate 









CONTROL REGISTER-TIMING GENERATORS CARD "C' 


F rom 

To 

Function 

Condi ti 

C3 

A3-G44 

1 us parallel load 

output 

C5 

D5 

serial shift gated clock 

output 

C49 

D49-10's tws pole 

sequence counter stop sweep 

; npul 

C51 

D54 

15 ms T.D. start sweep 

input 

C53 

D53-N.O. SW 

manual start pushbutton 

input 

C55 

N.C. SW 

manual start pushbutton 

input 

C57 

G30 

start sweep 

output 

C4 

D15-A5 

8p s shift clock 

output 

C6 

GI 6 

cryogenic cavity 

input 

C36 

D50 N.C. SW 

manual reset pushbutton 

input 

C40 

D28 

dummy data enable 

output 

C42 

D30 

real data enable 

output 

C52 

A52 

2qs shift register mode control 

output 

C58 

F55 

recorder clock 

input 


Cl 

7495 


C2 

7495 


C3 

7495 


C4 

848 


C5 

7404 


C6 

848 


C7 

RC network for one 

shot 

C8 

RC network for one 

shot 

Cll 

7400 


Cl 2 

848 


Cl 3 

9601 


C14 

9601 


Cl6 

7486 


C19 

832 


C26 

848 


C29 

7400 


C30 

848 



4-bit right/left shift register 
4-bit right/left shift register 
4-bit right/left shift register 
clocked flip-flop 
HEX inverter 
clocked flip-flop 
(Cl 3) 

(C14) 

quad 2-input positive NAND gate 

clocked flip-flop 

one shot multivibrator 

one shot multivibrator 

quad exclusive OR gate 

dual buffer 

clocked flip-flop 

quad 2-input positive NAND gate 

clocked flip-flop 








MODE COUNTER-SEQUENCE COUNTER- ANTENNA MULTIPLEXER 

CARD "C" 


F rom 

To 

C33 

10's tws pos 3 

C37 

10's tws pos 5 

C39 

10's tws pos 6 

C45 

10's tws pos 9 

C6 

D47 

C8 

N.O. sw 

CIO 

D10-G44 

C 1 2 

N.C. sw 

C 1 4 

N.O. sw 

C 1 6 

G46 

C 1 8 

G48 

C20 

G50 

C22 

G52 

C24 


C26 

A36 

C28 

A38 

C30 

D20-E20-G40 

C32 

DI6-E16-G36 

C34 

DI2-E12-G32 

C36 

D50-N.C. sw 

C44 

D48 


Function 

sequence counter 30 counts 
sequence counter 50 counts 
sequence counter 60 counts 
sequence counter 90 counts 
cryogenic cavity gated 
manual reset pushbutton 
reset and antenna multiplex 
antenna multiplex enable 
antenna multiplex enable 
mode 1 reset 
mode 2 reset 
mode 3 reset 
mode 4 reset 
mode 5 reset 
antenna identification 1 
antenna identification 2 
antenna driver 1 
antenna driver 2 
antenna driver 4 
manual reset pushbutton 
manual antenna advance 


Condition 

output 

output 

output 

output 

input 

input 

input 

input 

input 

output 

output 

output 

output 

output 

output 

output 

output 

output 

output 

input 

input 
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w f 


7490 

7400 

7400 

7400 

7486 

848 

848 

cambion 

7490 

7400 

7400 

846 

7490 


decade counter 
quad 2- input NAND gate 
quad 2-input NAND gate 
quad 2-input NAND gate 
quad exclusive OR gate 
clocked flip flop 
clocked flip flop 
for wire connections 

decade counter 

quad 2-input positive NAND gate 
quad 2-input positive NAND gate 
quad 2- input NAND gate 
decade counter 
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GENERAL PURPOSE TIMING EVENTS CARD "D" 


F rom 

To 

Function 

Condi ti on 

D3 

D5 

C5> 

5us serial shift gated clock 

input 

D15 

A5-C4 

8|js serial shift gated clock 

output 

D19 

F16 

end of file Sps 

output 

D41 

N.O. sw 

end of file pushbutton 

input 

D43 

N. C . sw 

end of file pushbutton 

output 

D47 

C6 

cryogenic cavity gated Sms 

output 

D6 

G16 

cryogenic cavity gated 

i nput 

DIO 

CIO 

reset, antenna multiplex and ISmsTD 

input 

D48 

A54-G54 

real time counter reset 

output 

D50 

C36-N.C. sw 

manual reset pushbutton 

1 nput 

D54 

C51 

1 Sms time delay 

output 
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D 1 

D2 

D3 

D5 

D6 

D7 

D8 

D9 

DIO 

Dll 

D12 

D18 

D24 


RC network for one shot (10uf, 10 K) 
RC network for one shot (.OOluf, 10 K) 
RC network for one shot ( . 00 5(_i f , 10 K) 
RC network for one shot (.OOluf, 10 K) 
RC network for one shot (10uf, 10 K) 


9601 

one shot multivibrator 

9601 

one shot multivibrator 

9601 

one shot multivibrator 

848 

clocked flip-flop 

9601 

one shot multivibrator 

9601 

one shot multivibrato 

7400 

quad 2- input positive ' AND gate 

7486 

quad 2- input exclusive OR gate 
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GENERAL PURPOSE FUNCTIONS CARD "D" 


* 

i 


F rom 

To 

Dll 

FI 1 

D13 

led 1 panel 

D17 

FI 5 

D35 

F35 

D37 

F55-C58 

D49 

C49 10's tws 

D51 

F53 

D53 

C53-N.O. sw 

D55 

led 2 panel 

D57 

N.O. sw 

D8 

C44 

D12 

El 2-G32- C34 

D 1 , 

led ant "0" 

D16 

E16-G36-C32 

D18 

led ant "1" 

n?.o 

E2n-G40-C30 

D22 

led ant "2" 

D24 

N.O. sw 

D2b 

N.C. sw 


Function 

Condition 

tape running 

input 

tape running ; ndicator 

output 

end of tape 

input 

write command 

output 

write shift clock 

input 

sequence counter switch pole 

input 

file protect missing 

input 

manual start pushbutton 

input 

file protect missing indicator 

output 

recorder on clamp 

input 

manual advance antenna 

output 

antenna "0" 

input 

antenna M 0" indicator 

output 

antenna "1" 

input 

antenna "1" indicator 

output 

antenna ”2" 

i nput 

antenna "2" indicator 

output 

manual antenna advance 

input 

manual antenna advance 

5 nput 


) 

} 




. % »*-« - , t * rr^VJT- , *v> . < ■ ' <■ - . » -,•*.*— * • MV'.* ^ ( 1V M-M " * * • -' .r . |K|£ m» ^r^-' 


D18 

7400 

quad 2-input positive NAND gate 

D19 

7403 

quad 2- input positive NAND gate 

D20 

848 

clocked flip-flop 

D24 

7486 

quad 2- input exclusive OR gate 

D25 

510Q resistors 


D26 

7403 

quad 2-input positive NAND gate 

D27 

7400 

quad 2-input positive NAND gate 

D28 

848 

clocked flip-flop 
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ANTENNA SWITCH DRIVERS CARD "E” 


8 

[ 

I 

S 


Function Condition 

+ 28 volts input 

antenna "0'' output 

antenna "0" input 

antenna "1" output 

antenna ''1" input 

antenna "2" output 

antenna "2" input 


i 











CARD "G" 


•# 


From 

To 

G6 

connector 1 

G8 


G10 


G14 


G16 

D6 

G18 


G20 

connector 2 

G22 


G24 

D52 

G54 

D10-C10 

G56 

D58-A58 

G58 




% 


Function 

Condition 

cryogenic cavity signal 

input 

amplifier "1" 

test point 

amplifier "2" 

test point 

amplifier "3" 

test point 

cryogenic signal gated 

output 

cryogenic comparator signal 

test point 

reference cavity signal 

input 

reference cavity amplifier "1" 

test point 

reference cavity comparator 

output 

mode counter reset signal 

input 

100 K Hz clock 

input 

reference cavity amplifier "2" 

test point 
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G1 

G2 

G4 

G5 

G7 

G8 

G9 

G10 

Gil 

G13 

G14 

G15 

G20 


LM301 operational amplifier 
1MQ pot and lMfl resistor 
LM301 operational amplifier 
lMfi pot and lMfi resistor 
LM301 operational amplifier 
capacitor cambion 

848 clocked flip-flop 

LM301 operational amplifier 

capacitor cambion 

LM301 operational amplifier 

LM301 operational amplifier 

7400 quad 2- input NAND gate 

848 clocked flip-flop 









F rom 

To 

Function 

Condi tion 

G28 

E4 

+ 28 volt 

output 

G32 

D12-E12-C34 

antenna ’O’’ 

i nput 

G34 

E10 

antenna M 0” drive 

input 

G36 

D16-E16-C32 

antenna "1" 

input 

G38 

E14 

antenna " 1 " drive 

input 

G40 

D20-E20- C30 

antenna ”2" 

input 

G42 

E18 

antenna ”2" drive 

input 

G44 

C3-A3 

parallel load 

input 

G46 

C16 

mode 1 counter 

input 

G48 

C18 

mode 2 counter 

input 

G50 

C20 

mode 3 counter 

input 

G52 

C22 

mode 4 counter 

input 

G54 

D48-A54 

15 ms time delay reset 

input 
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G3 

cambion for wire 

connections 

G6 

cable cambion 


G12 

7400 

quad 2-input positive NAND gate 

G17 

cambion for wire 

connections 

G18 

7404 

HEX inverter 












6 STAGE DECADE COUNTER CARD 


" 

2/28/73 Connector Pins 

Function 

f 

1 

+ 5 V 

* 

2 

+ 5 V 


3 

externa 1 reset 


4 

signal input 


5 

strobe input 

r 

6 

segment test 

4 

56 

push button reset (N. O. ) 

t 

57 



58 

push button reset (N. C. ) 


59 

common 

4 

60 

common 


! 


5 , 

fe 

i 


2/28 / 73 


INTERFACE CARD 


Connector Pins 

Function 

1 

+ 5 V 

2 

+ 5 V 

3 

100 KHz clock 

4 

common 

5 


6 


7 


8 


9 

strobe #1 ANT "0" 

10 

common 

11 

strobe #2 ANT "1" 

12 

common 

13 

strobe #3 ANT "Z" 

14 

common 

15 

external reset (all resets 

16 

common 

17 


18 


19 

parallel load 

20 

common 

59 

common 

60 

common 
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Table E-I. Interface Connector Pin Li-t 


■ 

F 


CTIVE 

WIRE 

OLOR 

\CT. 

S1G. 

SIGNAL (see notes 1 through 3) 

PIN 

GRD. 

1— 


< 1 



<7> 


TRANSPORT CONTROL SIGNALS 


3 


Run Normal Speed 

input level 

80 

X 

X 

X 

X 

V 

rr 

' * 

X 

n 

5 


Run High Speed 

input level 

30 

X 

X 

X 

X 

A 

X 

X 

j 

7 


Foi ward Select 

Input level 

70 

X 

X 

X 

X 

X 

X 

X 

1 

9 


Reverse Select 

input level 

6C 

X 

X 

X 

X 

X 

X 

X 


1 1 


Tape Running 

output level 

£0 

X 

X 

X 

X 

X 

X 

A 


13 


Load Point 

output level 

8 

X 

X 

X 

X 

X 

X 

X 


15 


End of Tape 

output level 

6 

X 

X 

X 

X 

X 

X 

X 


17 


Broken Tape 

output level 

5 

X 

X 

X 

X 

X 

X 

A 


19 


Unload Command 

input pulse 

10 

X 

X 

X 

X 

X 

X 

X 


21 


Rewind to Load Point 

input pulse 

90 

X 

X 

X 

X 

X 

X 

X 


23 


Off Line 

input pulse 

9 

X 

X 

X 

X 

X 

X 

X 


25 


Rewind in Process 

output level 

2d 

X 

X 

X 

X 

A 

X 

X 


27 


On Line 

output level 

40 

X 

X 

X 

X 

X 

X 

k 




Reserved - Oo Not Use 


- 








X 



Reserved - Do Not Use 


0 








X 

29 


Chassis Ground (see note 4) 


1 

X 

X 

X 

X 

v 

JL 

X 

X 


Wdrrc- 


iLS 








31 


Write Select 

input level 

7 

X 

X 

X 




X 


33 


Write Parity Select (odd) 

Input level 

9? 

;< 



X 


X 

X 




Reserved 


0 

X 

X 






X 

35 


Write Command 

Input pulse 

50 

X 

X 


X 


X 

X 


18 


Data Channel 0 (Write) 

input level 

97 


X 


X 


X 

X 


20 


Oats Channel 1 (Write) 

Input level 

90 


X 


X 


X 

X 


37 


Oats Channel 2-B (Write) 

input level 

93 

X 

X 


X 


X 

X 


39 


Data Channel 3-A (Write) 

input level 

96 

X 

X 


X 


X 

X 


41 


Data Channel 4-8 (Write) 

Input level 

95 

X 

X 


X 


X 

X 


43 


Data Channel $-4 (Write) 

input level 

94 

X 

X 


X 


X 

X 


45 


Data Channel 6-2 (Write) 

input level 

8 

X 

X 


X 


X 

X 


47 


Data Channel 7-1 (Write) 

Input level 

98 

X 

X 


X 


X 

X 


49 


Data Channel P-C (Write) 

input level 

30 

X 

X 


\ 




X 

51 


Write Status 

output level 

90 

y 

X 


X 


X 

X 


S3 


File Protect Ring Missing 

output level 

98 

X 

X 


X 


X 

X 


55 


Write Clock 

output pulse 

70 

X 

X 


X 


X 

X 


57 


Write Clock Gate 

output level 

SO 

X 

X 


X 


X 

X 


22 


Write Echo Error 

output pulse 

10 

X 

X 


X 


X 



12 


2X Write Clock 

Input pulse 

4 

X 

X 






X 

14 


Write ETC 

Input pulse 

20 

X 

9 


X 


X 

X 


16 


EOF Command 

input pulse 

80 

X 

X 


X 


X 

X 
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APPENDIX C 


Q MEASUREMENT SUMMARY 

The objective of the Q measurement task was to determine che 
spherical vessel Q as a function cf diameter and wall material. The 
results of this measurement task are summarized in tabular form. 


Spherical Vessels 


Mode 

19” Cu 

18" SS 

60" SS 

™,i 

91,000 

1 1, 600 

6, 200 

™21 

85, 200 

10, 900 

13,600 

TE u 

81, 100 

19, 700 

28,400 

TM, 

< 1 

41, 200 

10, 900 

17, 500 


Considering the lowest Q (TM^ in the 60” SS), f^ is 172 MHz when 
empty and would be 155 MHz when filled with boiling point LH^. The 
total frequency shift from empty to full would be 1 7 MHz, The band- 
width at the 3dB points is 0.028 MHz, or approximately 0. 15% of the 
empty-full bandwidth, For LOX the 3dB bandwidth is approximately 
0. 10% of the empty-full bandwidth. Thus it appears that the Q will not 
seriously degrade the RF resonance technique in stainless steel spheri- 
cal vessels up to 5 feet in diameter. 

Let us assume that the change in Q pattern measured with the i8 n 
and 60 M stainless steel vessels continues to even larger diameters. 

Then for a 15 foot diameter stainless steel vessel, the 3dB bandwidth 
is 0.3% of the empty full LH 2 bandwidtn. The LOX 3dB bandwidth is 
approximately 0. ?% of the empty -full bandwidth. 
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Figure C3. Frequency response of the copper spherical vessel for the 
TEu mode. 
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Figure C4. Frequency response of the copper spherical vessel for the 
TM 031 mode. 
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Figure C5. Frequency response of the 60 inch diameter stainless steel 
sphere - TMqh mode. 
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] Figure C6, Frequency response of the 60 inch diameter stainless steel 

' sphere - TE 0 n mode. 
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Figure C7. Frequency response of the 60 inch diameter stainless steel 



APPENDIX D 


DATA REDUCTION FROM THE MAGNETIC TAPE UNIT 

The magnetic tape unit stores the information coming from the 
data acquisition system described in Appendix B. Each data run lasts 
about ten seconds; during this time the system records about 25 inde- 
pendent measurements for each mode or about 375 data points in all 
(counting 3 antennas and 5 modes for each antenna). The data reduction 
processes the raw tape, finds the average and standard deviation of the 
25 data points for each mode and punches these numbers on to standard 
data cards along with other identification codes. Data from other 
sources are also punched on these cards i. e. , mass, pressure and 
temperatures. The cards are then ready for plotting routines, and 
other types of data analysis. 

I Program DMPMODES 

Since the tape recorder in the laboratory is a continuous write 
recorder and not a sequential recorder, the letter D is hard-wired to 
the output whenever there is no data being recorded. Most of the tape 
is, in fact, dummy data. DMPMODES accomplishes many functions. 

It takes the data from the magnetic tape made in the laboratory, deletes 
all D 1 8 and other non-essential characters, sorts the 12 character data 
by antenna, lists this data and writes this data on a second magnetic tape. 

Essentially, the program works as follows. A tape record is 
read into array M. A pointer, KNT, moves down the array. If the 
word in M(KNT) is all D's, KNT is incremented. When a word that is 
not all D's comes along, the character: (or 15 octal) is searched for. 
When a: is found, it and the characters surrounding it are placed into 
the first two locations of the array IW. Thinking of these two 8 byte 
words as one 16 byte word, it should look like < ^i^2 r l r 2 a:t l t 2 t 3 t 4 t 5 t 6 
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DDDD, where all subscripted letters are integers, and D's are 

dummy data, d^d^ * s t ^ e data code selected from the recorder front 

panel. r \ r 2 * s run cot ^ e > a l so selected from the front panel. 

a . is the antenna code, t -t, are the data. If this 16 byte "word" is 
I 1 o 

not correct, that is, if an alpha character appears where an. integer is 
supposed to be, the data point is rejected and the program starts search- 
ing for the next data point. If the data point appears to be all right, it is 
decoded, and the data is printed under the appropriate antenna column. 
Also, the 16 byte "word" is written on a second magnetic tape, thus 
saving later programs from having to re-sort the mass of characters 
on the input tape. 

II. Program AVMODES 

Basically, this program takes data from the second magnetic 
tape, sorts it by antenna and mode, finds the average time and standard 
deviation for each mode for each antenna, lists these averages and 
standard deviations, and punches on standard data cards the averages 
and standard deviations, along with the date, date code, run code, 
antenna number, mode number, and the number of points used to calcu- 
late the average and standard deviation. The input magentic tape is 
essentially a sequence of the magnetic tapes written by DMPMODES. 

The tape has all the data from all the days of running the experiment. 
Each day is separated from the others by an end of file marker, and 
each run is separated from other runs by an inter-record gap. The 
first record of each file contains the date of the run. 

The program works as follows. The first record is read into 
the first 10 words of array M. The date of the run is extracted from 
these ten words and stored into the variable IDATE. Each of the other 
records in the file are read into the array M, one at a time. Every 
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two words of the array M are the same as in the 16 byte "word" described 
above in the discussion of DMPMODES. This 16 byte word ; <? decoded, 
or broken up, into the date code, IDC ; the run code, I RC; ti.c antenna 
code, LA.; the colon character, ICOL; and the time data, ITIME. The 
time is checked to see which mode, if any, it belongs to. The point is 
then used to calculate the average and standard deviation of that 
particular mode of that particular antenna. Note that for each record 
on the tape, or run, there is a possible total of 1 5 of these averages 
and standard deviations to be printed and punched (corresponding to three 
antennas with five modes eaci'l. The punch card output is used for 
graph routines and other data analysis programs. 

III. The following pages give a listing of these programs which are 

being used on the CDC 3800 computer. 
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program dmpmooes 

OIMENS I ON MSG ( 3 ) » I W ( 8 I * HZ I J ) * I HZ SO ( 3 ) . ' < 3 ) . T SO( 3 ) .X ( 3 ) . IC*. 5000) 
COMMON M( 32765) 

1 FORMAT ( 15 ) 

2 FORMAT (*1REC0RD *.I5»* LONGER THAN 32764* SOME LOST.*) 

3 FORMAT (8R1) 

4 FORMAT (2!2*IltRltl6.4X) 

5 FORMAT <*1*.57X.*TAPE RECORD NUMBER *,I3.3A8/*8* 

I51X«*DATE CODER * . I 2 ♦ 5 X **PUN CODER **I H 

255X **ANTFNNA IDENTIFICATION CODE* / 3 3X * *0* ♦ 32 X ** 1 * ♦ 32 X .*2* > ) 

6 FORMAT ( 3 1 X * 16 • 

7 FORMAT I 64X « 1 6 ) 

8 FORMAT ( R7X * 1 6 ) 

9 FORMAT (*CAVG*»2 7X,F9.2*2I23X*F9.2 •/* S IGMA* . 23X . E 1 2. 6 * 2 < 2 1 X .£ 12 .5 > 
1/1X.I4,* SCANS DELETED IN THIS RECORD.*' 

NRUNO*RRR 
CALL 09OVFR 
CALL iohcheck 
N9EC *0 

READ 1.NSKIP 
IFICOF. 6010. 103 
103 DO 105 I * 1 *NSK I P 
BUFFER JN< 1.0) (M.M) 

IQ* IFIUNIT.l ) 104,105 
105 CONTINUE 

10 BUFFERIN< 1,0) (M.MI32765) ) 

IFINREC.EO.O.OP.XCA.EQ.llGO TO 115 
BUFFEROUTl 2*1 > ( IOA . IOA ( rOA-l ) ) 

115 KOA-1 

MSG 1 1 >*MSG<2 ) -MSGI 3 )*8H 
KNT “0 

INZ C 1 ) » IHZ ( 2 • ■ IHZ ( 3 ) *0 

K « 1 « 2 I «K « 3 » *C 

IHZSO ( 1 ) ■ ’ HZ SO ( 2 ) « I HZSO ( 3 ) «0 

NREC-NRFC*2 

IOROP— 1 

LINE-75 

KCODE-O 

117 IFlUNlT.2lli7.il 

11 IFIUNIT.l 111.15.99.13 
13 MSG ( 1 >-8H CONTAIN 

MSG( 2 ) »8H*i PARITY 
MSG( 3 ) -8H ERRORS 
15 L-LENG't HF 1 1 ) 

IFIL.GT.4 )f»0 T 0 195 
NREC-NREC-1 
GO TO 10 

193 IF I L .L T .32764 )GO TO 198 
PRINT 2.MREC 
198 l DROP- I DROP* 1 
203 KNT-KNT* 1 


IFUCODE.EO.l KnT.KMT + 1 
KCODF *0 

20 IFUNT.GT.L-2 »23*24 

23 DO 242 1=1*3 
F I »K C I ) 

T ( I 1*1 HZ < I > 

T ( I ) «T ( I ) /F I 

IFUt P .GT.l )GO TO ?41 

TSO( I > *0 • 

GO TO 242 

241 TS0( I l = IHZSO < I ) 

TSQ( I )=SQRT(ABS( (FI*TSO( 1 )-T « !) *»2 > /F 1 / < F I -1 . ) )) 

24 2 CONTINUE 

PRINT 9*1 T( I ) * 1*1.3 ) * < TSO « I > » 1=1 *3* . I DROP 
GO TO 10 

24 IF(MOCNT) .EQ.8HDDDDDDDD)GO TO 27 
DECODE (8*3*M(KNT) MW 

25 DO 26 1=1*8 

IF* IW< P.NE.15B>GC TO 26 
NSH« I -6 

I F ( NSH >29*37*30 

26 CONTINUE 

27 KNT=KNT+1 
GO TO 20 

29 knt =kn r« 1 

NSH*N5H+8 
KCODE« 1 

30 MSK1«2**< 6*NSH) 

MMNl= M SK 1-1 
M SK2»?** I I 8-NSH ) *6 1 
MMN2=MSK2~1 

DO 35 1=1.2 

IW( p«MUNT+I-l ).AND.VMN2 
IW(P = IW< I»*MSK1 
IC-MIKNT+I I/MSK2 
35 IW< I> = IW< P»OR. ( IC.AND.MMN1 • 

GO TO 40 
37 !W<P=M(KNT> 

I W I 2 >«m{(cnT + 1 ) 

40 ITEST«1W(2UAND.770000000QB 

IF ( ITEST •EO.0KGOODOOOO >G0 TO 198 
ITEST= IWI p .AND. 77000000000000009 

IFI ITEST.EQ.8HD0000000) IW( 1 i « IW < 1 ) .AND. 777777777777770 
DECODE 1 16*4* I Wl ID.NRUN.NANT* 1C* I D^TA 
IF«' I OHERRIOP 198*43 
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43 I F ( IC.NE. 158)60 To 1 9P 
IOA(KOA). I W C 1 ) 

IOAUOA + l >'IW(2> 

K0A*K0A+2 
LINE-L INE+1 

IF ( NRUN.EC.NRUNO.AND.L INE.LE .67 >60 TO 53 
PRINT 5»NREC*(MSG( I ) • 1 = 1*5) » I D * NRUN 
LINE=1 

53 NANT*NANT+1 

K»NANT )*K (nANT)+1 
IHZ(NANT)sIHZ(NANT)* I DAT a 
IHZS0(NANT)*IHZ?.C(NANT )♦ I DATA* I DATA 
NRUNG»NRUN 

GO TO (55*56*57) NANT 
*5 PRINT 6.IDATA 
GO TO 205 

56 PRINT 7 ♦ I DATA 
GO TO 205 

57 PRINT 8 ♦ I DAT A 
GO TO 205 

v« STOP 
END 



PROGRAM avmodes 

0 IMENS ION M ( 5000 I , T I ME < 5 *3 > * T I MESQ ( 5 . 3 ) .KNT ( 5 * 3 > .MODE < 5 . 2 > 

DATA! ( (MODE! I.J»»I*1*5J.J=1*2 1=1100.2900.4200.5100.6500.1655.3900. 
15100.6060.7800) 
cn _ IOHCHECK 
NF : 0 

RE A r ’ 1 .NFSKlP 

1 FCkMA'i t 15 ) * 

I F { EOF .60 ) 10.10? 

102 DO 105 1 = 1 «NF SKIP 
105 CALL SKIPFILEU ) 

10 BUFFER INt 1.1 > <M,M( 10) > 

NREC=0 

NF = Nr + 1 

11 IF < UN I T . 1 ) 11.15.99.13 
13 PRINT 2.NREC.NF 

2 FORMAT (*0PARITY ERROR IN RECORD *»I3* # OF FILE *.I3> 

15 IDATE= M < 1 ) .AND. 7777777777B 
IDATE= IDATE*1000000B 

M ( 2 )=M ( 2 ) / 1000000000 OB 
I DAT E= I DATE. OR. (M( 2 ) .AND.777777B ) 

PRINT 6 

6 FORMAT ( 1H1 / ) 

16 DO 1’ J= 1 *3 
DO 17 1=1.5 
KNT( I.j>=0 
TIME( I » J ) =0 • 

17 T IMESQ ( I . J ) =0 • 

BUFFERINI 1.1 ) t M »M( 5000 ) ) 

NREC=NREC+1 

21 IFIUNTT.l >21.25.10.23 
?3 PRINT 2.NREC.NF 

25 L«LENGTHF(1> 

DO AO J=1 *L*2 

3 FORMAT (2I2.I1.R1.I6.4X) 

DECODE! 16. 3.M(J>) I DC. IRC. I A. I COL. I TIME 
IF ( IOHERR ( 0 * ) AO .26 

26 IF ( ICOL.NE. 15B >G0 TO AO 

IF! IA.LT.O.OR.IA.GT.2>GO TO AO 

IAC-IA+1 

DO 30 1=1.5 

IF ( IT I ME. GT. MODE! I .2 ) ) GO TO 30 
f c ( I T I ME . GE • MODE ( I ♦ 1 • ) GO TO 32 
GO TO AO 
30 CONTINUE 
GO TO AO 

32 KNT < I » I AC ) «KNT ! I ♦ I AC • +1 * « 

PN*ITIME • 

PN-PN/100. 

TIME! I . I AC ) =T IME ! I ♦ I AC ) *PN , 

T I MESQ! I » I AC ) «T IMESQ! I . IAC > +PN»*2 
40 CONTINUE 

50 DO 55 J«l»3 ] 

I A» J- 1 ^ 

4 


1 
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00 55 1*1*5 

PN=KNT ( I , J) 

IF ( <NT ( I*J».EfJ.0>G0 TO 55 
IFIKNT ( I * J) .GT.i )r-0 TO 52 
T IMESO < I»J>=0. 

GO TO 53 

52 T IMESO < 1 . J) =SORTl ABS IPN*T I MESQl I « J 1 -T I ME < I , J ) **2 ) / (PN*<PN~1 . ) ) ) 

53 TIMEC l «J)=TIME< I.JI/PN 

PUNCH 4.IDATE.IDC.IRC. IA»I .TIME! I . J> *TIMESG( I *J) *KNT( I , J) 

4 FORMAT < A8.213.2 12 . F? .3 * E12 .4 * 1 5 > 

PRINT 5»IDATE*!DC* IRC. I A* I *TJMF.< I *J> *TIMESQ( I *J> *KNT( I * J ) 

5 FORMAT (1X*A8*2I3.2I2»F7.3*E12.4.I5> 

55 CONTINUE 

GO TO 16 
99 STOP 
END 
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APPENDIX E 


UNIFORM DENSITY HYDROGEN 
DATA ANALYSIS AND ACCURACY STATEMENTS 

The data analysis for uniform density hydrogen is based on 41 
observations at tank pressures above 400 psi; the temperatures at the 
top of the tank were slightly lower than the bottom temperatures (the 
heater was at the bottom) indicating almost uniform temperature with a 
condition for convective mixing due to gravity; experimentally, these 
conditions were necessary to achieve a near uniform density within the 
experimental vessel. 

The data for the TM^^ mode, as plotted in Figure 6 of the main 
text, is expected to follow equation (6), i. e. , 



(E- 1 ) 


Using equation (2), the expression for f^ 


in the time domain 


f l = f 0 + r(t l " V 


(E-2) 


it follows after a little algebra 



(E-3) 


or equivalently 


E-l 


JL. 


J * 


where 


* A. 01 

At. = At + 

1 Or 



2AM 


“l = *1 - ‘o 


at = f oi ' f o 

r 


(E-4) 


(E-5) 


(Atg is just the time interval at M = 0). The gauging function is then 
between At and M; these quantities being related to three independent 

1 fo i a 

parameters At , , — through equation (E-4) . These three param- 

Or V 

eters may be determined by fitting the data to (E-4) or, alternately, 
may be determined by separate physical measurement. For example, 
the sweep rate, r, can be determined from Atg, fg and fgj by 

f rtl - f 
01 0 


The measured values of these quantities were fg^ = 581.9 MHz, 

fg =411 MHz and Atg = 16.20 milliseconds giving r = 10,54 MHz/msec. 

V can be calculated from the inside radius of the tank R = 8. 906 

4 3 

inches (22. 62 cm) giving V = 4. 85 x 10 cm . A can be obtained 

from other experimental data; ref. 4 of the main text gives A = 1 . 006 
3 

cm / gm. Table E-I gives the comparison between these values and the 
values obtained by a nonlinear least squares fit of the data to equation 
(E-4).* 


* C. Daniel, F. S. Wood. Fitting Equations to Data, Wiley Interscience, 
N. Y. (1971) p. 320. 
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TABLE E-I 




r 


A 

V 


Measured value 16.20 msec 


55.15 msec 


9. 42 x 10 



1 


I itted value 16.19 


52. 88 


9. 44 x 10 


95% confidence 

limit (lower) 16.20 39.20 6.81 x 10 


95% confidence 

limit (upper) 16.20 66.50 12.21 x 10 


Although the confidence limits for the fitted values are relatively 
large, the actual fitted values are quite consistent with the independently 
determined values. 


The semi -empirical gauge equation may then be given by 
At, = 16. 19 + 52. 88 


l - 9. 44 x 10 M 


1 


1 + 2 x 9. 44 x 10 -3 M 


■ 


(E-7) 


where At^ is measured in milliseconds and M in pounds; it remains 
to determine the accuracy of the data with respect to this equation. 

Table E-II gives the observed value of At^ and the value of At^ 
calculated from (E-7) for each of the observed mass values. The 
residuals are plotted in Figure E-l as the Percent Full Scale Residual 
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Mass (lbs) At (OBS) 

0.00 16.200 

0.20 16.018 

0.44 15.839 

0.52 15.820 

0.73 15.631 

0.84 15.550 

0.94 15.460 

1.00 15.500 

1.15 15.340 

1.15 15.329 

1.16 15.326 

1.42 15.163 

1.42 15.160 

1.54 15.034 

1.54 15.029 

1.56 15.029 

1.57 15.019 

1.58 15.030 

1.81 14.842 

2.06 14.674 

1.85 14.841 

2.11 14.670 

2.43 14.415 

2.45 14.413 

2.82 14.180 

3.33 13.820 

3.70 13.557 

4.23 13.141 

4.25 13.136 

4.71 12.815 

5.19 12.501 

5.20 12.518 

5.62 12.221 

5.63 12.222 

5.68 12.212 

6.43 11.679 

6.44 11.686 

6.51 11.680 

6.77 11.480 

6.79 11.480 

7.08 11.304 


TABLE E-II 
At ^ (CALC) Residual 

16.186 0.014 

16.037 - 0.019 

15.859 - 0.020 

15.799 0.021 

15.644 - 0.013 

15.564 - 0.014 

15.483 - 0.023 

15.446 0.054 

15.337 0.003 

15.337 - 0.008 

15.330 - 0.004 

15.141 0.022 

15.141 0.019 

15.054 - 0.020 

15.054 - 0.025 

15.039 - 0.010 

15.032 - 0.013 

15.025 0.005 

14.859 - 0.017 

14.687 - 0.013 

14.831 0.010 

14.645 0.025 

14.417 - 0.002 

14.403 0.010 

14.142 0.038 

13.786 0.034 

13.530 0.027 

13.167 - 0.026 

13.154 - 0.018 

12.856 - 0.041 

12.520 - 0.019 

12.514 0.004 

12.235 - 0.014 

12.228 - 0.006 

12.195 0.017 

11.703 - 0.024 

11.697 - 0.011 

11.651 0,029 

11.483 - 0.003 

11.470 0.010 

11.283 0.021 


% Full Scale Residual 

0 . 28 % 

- 0 . 38 
- 0.41 

0 . 43 
- 0 . 26 
- 0 . 29 
- 0 . 47 

1 . 10 
0 . 06 

- 0 . 16 
- 0 . 08 
0 . 45 
0 . 38 
- 0 . 41 
- 0 . 51 
- 0 . 20 
- 0 . 26 
0 . 10 
- 0 . 35 
- 0 . 26 
0 . 20 
0 . 51 
- 0 . 04 
0 . 20 
0 . 78 
0 . 69 
0 . 55 
- 0 . 53 
- 0 . 38 
- 0 . 83 
- 0 . 39 
0 . 08 
- 0 . 28 
- 0 . 12 
0 . 35 
- 0.49 
- 0 . 22 
0 . 59 
- 0 . 06 
0. 20 
0.43 
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Figure El. Residuals vs At 




vs. Atj(calc); this graph shows a fairly random scatter of the data 
between -0. 83% and 1. 10% full scale. 

The cumulative distribution of the residuals is plotted in Figure 
E-2 vs. the expected intervals in a normal distribution. The S-shaped 
trend of the data indicates that the frequency distribution is slightly 
flatter than the normal bell- shaped error function. This is probably 
due to the fact that the At^ are really only measured to four significant 
figures, the most significant digit corresponding to 0. 2%; this would 
cause a broadening of the frequency distribution on the order of ■* 0. 2%. 
The slope of the straight line through the tails of the curve gives an 
estimate of the standard deviation a = 0.4%. This gives a 3a deviation 
(the 99.9% confidence interval) of 1. 20%. 

Operational Readout 

As indicated in the main text, a quadratic fit will be easier to 
work with operationally and should give a sufficiently accurate gauging 
function. The quadratic equation 

At = 16.20 + AM + BM 2 (E-8) 

was fitted to the data in Table II by solving (E-8) simultaneously for 

M = 6.77, At = 11.48 and M = 2.43, At = 14.415. The result is 
1 _ 3 1 

A = -0.7554 and B = 8.61 x 10 . The residuals for this fit range 

between -0.63 and 0.93 percent full scale, a slightly better fit than 
the theoretical curve. This analysis is tabulated in Table III. 


* Chernov, H. and Lieberman, G. S. , J. Amer. Statist. Assn. 49, 
778-85 (1954). 
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TABLE E - III 


Mass (lbs) 

A tj ( OBS ) 

At i (CA I 

0.00 

16.20 

16.2 

0.2 

16.018 

16.049 

0.44 

15.839 

15.869 

0.52 

15.82 

15.809 

0.73 

15.631 

15.653 

0.84 

15.55 

15.571 

0.94 

15.46 

15.497 

1.00 

15.50 

15.453 

1.15 

15.340 

15.342 

i .15 

15.329 

15.342 

1.16 

15.326 

15.335 

1.42 

15.163 

15.144 

1.42 

15.160 

15.144 

1.54 

15.034 

15.057 

1.54 

15.029 

15.057 

1.56 

15.029 

15.042 

1.57 

15.019 

15.035 

1.58 

15.03 

15.027 

1.82 

14.842 

14.853 

2.06 

14.674 

14.680 

2.11 

14.67 

14. 644 

2.43 

14.415 

14.41 5 

2.45 

14.413 

14.400 

2.82 

14.180 

14.138 

3.33 

13.820 

13.780 

3.7 

13.557 

13.522 

4.23 

13.141 

13.158 

4.25 

13.136 

13.145 

5.19 

12.501 

12.511 

4.71 

12.815 

12.833 

5.2 

12.518 

12.504 

5.62 

12.221 

12.226 

5.63 

12.222 

12.220 

5.68 

12.212 

12.187 

6.43 

11.679 

11.698 

6 . 44 

11.688 

11.692 

6.51 

11.68 

11.647 

6.77 

11.480 

11.480 

6.79 

11.480 

11.467 

7.08 

11.304 

11.283 


C) "'a Full Scale- Residual 

n. no'’:, 

0. e, } 

0.62 
0.21 
-0.45 
-0.44 
-0.77 
0.93 
- 0 . 0 ". 

- 0.30 
- 0 . 10 
0.37 
0.31 
- 0.47 
- 0 . 57 
- 0 . 28 
- 0.33 
0.04 
- 0.24 
- 0 . 13 
0.52 
0.01 
0.25 
0.85 
0.82 
0.70 
- 0.36 
- 0 . 19 
- 0.21 
- 0.37 
0.27 
- 0.12 
0.04 
0.51 
- 0.40 
- 0.09 
0.67 
- 0.01 
0.25 
0.42 


E-8 





APPENDIX F 


TOTAL MASS GAUGING IN A SPHERICAL RESONANT CAVjTY 
Int roduction 

When a closed metal container is excited by an RF antenna probe 
inserted through a hole in the container, theoretically there are an infinite 
number of excitation frequencies for which the container is strongly 
coupled to the antenna; this me^ns that energy can flow more freely 
between the antenna and the container at these resonant frequencies. 

The resonant frequencies correspond to standing wave patterns in the 
cavity which are called resonant modes. The wave pattern of the mode 
which occurs at the lowest possible resonant frequency is called the 
fundamental mode. This mode and the modes of the next few higher fre- 
quencies are called lower order modes. 

When the cavity is uniformly filled with a fluid, the resonant fre- 
quency changes because the velocity of propagation of the resonant 
standing wave, c = 1 //pe, depends on the dielectric constant, e, and the 
magnetic permeability, p, of the fluid. For example, in a spherical 

resonant cavity uniformly filled, the resonant frequencies, f , are given 

np 

by 

f = y ” P >— (F-l) 

np 2Trb/ge 

where b is the radius of the sphere, and n and p are subscripts which label 

the different modes (these will be explained in detail). The u are 

np 

eigenvalues of the modes and are obtained in the process of finding solutions 
to Maxwells equations. The resonant frequencies can then be related to 
total mass by using the Clausius -Mossotti relation 

pp = fti < F - 2 > 


F-l 



where P is the polarizability of the fluid which is a slowly varyinu func- 
tion of the fluid density, P, 

If the cavity is uniformly filled with a liquid the eigenvalues, u 

n P 

are just numbers independent of the fluid within the cavity and there is 
therefore a simple relationship between resonant frequency and total mass. 
However, if the container is only partially filled with liquid, the rest of 
the cavity being a vacuum or a gas, then the values of u will depend on 
a and € of the liquid, the and of the gas, and the geometry which 
the liquid takes within the cavity; the resonant frequency, then, is no 
longer an unambiguous function of mass but depends on the liquid geometry 
as well. This is because the standing wave patterns are distorted because 
of the boundary conditions at the liquid-gas interface. However, the 
resonant frequency of each partially filled mode does lie between the com- 
pletely empty and completely full values 

u u 

SR— . f < 2L ( f . 3) 

2nb/^It np Znb 

o o 

and varies continuously between these values as the cavity is filled. This 
suggests that the resonant frequency at least approximately indicates total 
mass independent of geometry. 

The purpose of this note is to investigate the geometry effects for 
a spherical cavity with spherical symmetry of the liquid gas interface. 

This geometry is similar to a "aero-g" formation with the liquid clinging 
to the walls and a gas bubble in the middle of the cavity. Tl ' reason for 
choosing this geometry is that it is one of the few examples of a partially 
filled cavity for which the Maxwell Equations can be solved in closed form. 
Even though this gometry is particularly simple, it does give a reasonable 


F -2 



indication of the uncertainty which may be involved when the geometry is 
not known. Numerical examples are calculated for cases in which the 
liquid is hydrogen or nitrogen. 

From a practical point of view, the spherical cavity is an ideal 
container geometry for this method of mass gauging. The reason for this 
is that the spherical symmetry of the cavity wall creates a degeneracy 
in the modes. That is. there are a number of standing wave patterns 
which have the same resonant frequency. This results in the fact that 
the distinct resonant frequencies of the lower order modes are wid. ly 
separated and minimizes the effect of mode crossing in a partially filled 
cavity. Mode crossing occurs when, for a particular liquid geomet./, 
the resonant frequency of a higher mode falls below that of a lower mode. 
For example, if the liquid is nitrogen, mode crossing between the first 
two modes is impossible and for the next few higher modes is quite 
unlikely; this is established from the table of eigenvalues, Table 1 on 
page F-12 and the inequalities expressed in (F-3). 

The relative independence of the lower order modes suggests that 
they can each be monitored independently. Since each mode has its own 
geometry in the standing wave pattern, it seems reasonable that the 
modes themselves may be used to at least partially determine the fluid 
geometry. (Mathematically the problem reduces to this: Given some of 
the eigenvalues of a boundary value problem, how closely can the eigen- 
functions be approximated. ) In fact, it will be shown that for the spherical 
symmetry considered in this analysis, that for a liquid of unknown density, 
both the location of the liquid-gas interface and the density (hence the 
total mass) can be determined uniquely if and only if five modes are 
monitored simultaneously. The reason for this is that each mode deter- 
mines exactly one independent relation between the resonant frequency 



of that mode and the fi\ . unknown parameters e, p, e > u , and a, 

r o o 

where r = a is the liquid-gas interface. For most applications it is suf- 
ficient to assume that e ^ & pss 1, leaving only two unknowns, namely 

o o 

e and a. In this case, two modes will uniquely determine the total mass. 

Solutions for Maxwells Equations in Spherical Coordinates 

When the cavity is resonating at an angular frequency, jl , the 
time phase of the electromagnetic field is the same at all points within 
the cavity. Hence, for a loss free cavity the electric and magnetic fields 
can be written as the real parts of Ee iajt and He 1 ^, respectively, where 
E ana H are vectors which depend only on the spacial coordinates. The 
source free Maxwell Equations can then be written 

curl E = - iuupH 
curl H - iuueE 
div eE = 0 

div pH = 0. (F-4) 

It should be emphasized at this point that only two assumptions have 
been made, the cavity is loss free and it is source free; in practice 
these are usually very good assumptions for calculating resonant fre- 
quencies. A third assumption which we will now make, may have to be 
justified more carefully in any given situation: we assume that there 
are two regions within the cavity, each of which have uniform density. 

The technical advantage of this assumption is that derivatives of p and 
e are not involved; the equation ( F-4) can be solved in each region where 
p and e are constant and the boundary conditions are then modified to 
include the liquid-gas interface. The boundary conditions can be written 

jeE*n, pH*n, Exn and Hxn continuous at each boundary point j (F-5) 
where n is the unit normal vector to the surface at that point. Since 
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div E = 0 and div H - 0, both E and H can be expressed in terms of 
vector potentials G and F, 

E - curl F 

H = curl G ( F-6) 

where the Maxwell Equations impose consistency conditions between 
G and F. Two independent solutions may be obtained by choosing a 
coordinate direction, say r, the unit vector in the radial direction and 
finding fields which are perpendicular to r. If E is perpendicular to 
r we say we have a TE (transverse electric) mode. This situation may 
be assured if F is chosen to be 

F = f? ( F- 7) 

where f is a scalar function of the spatial coordinates. In this case we 
have from (F-4) and (F-7) 

E = curl f? 

H = - ~~ curl curl fr . (F-8) 

t«u p. 

If H is perpendicular to r we say we have a TM (transverse magnetic) 
mode. This situation may be assured if G is chosen to be 

G = gr (F-9) 

where g is a scalar function of the spacial coordinates. In this case 
we have from (F-4) and ( F-9) 

E = curl curl g? 

uu e 

H = curl gr . (F-10) 

The general solution for E and H may be obtained by a superposition 
of (F-8) and (F-10) 
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curl f? + - curl curl g? 

tx e 


H curl gr - curl curl f? . 

ix p, 


(F - 1 1 ) 


To find equations which f and g satisfy, we consider the TE and TM modes 
separately. For the TM mode (F-10) and (F-4) imply that 


or 


curl E - - ix p. curl gr 


0 


curl (E 4- ix ugr) 

This last relation is satisfied only if 

E + ix p.g? = grad cp 


( F- 12) 


( F- 1 3 > 


for some scalar function cp. Substituting ( F-13) into the second of equa- 
tion (F-4) we have 


curl curl gr = x 2 pc gr + ixs gradcp. 
Using the vector equation 

curl curl gr = v 2 gr - grad (v-gr) 

and 

x 2 pe = k 2 

we find that g satisfies the following equations 

(V s + k 2 ) gr = 0 

V-gr = - ix e Cp . 


(F-14) 


(F-15) 


(E- 16) 


A similar argument for the TE mode shows that f satisfies the following 
equation 


(V a + k a ) f? = 0 

v-f$ = iu> ixv 4r. 


(F-17) 
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Equations (A- 16) and (A- 17) are equivalent to the scalar Helmholtz equa- 
tions for g/r and f/ r with standard solutions given by' : ' 

B (kr) L m (0, cp) , ( F- 18) 

n n 

where the L m (0, cp) are spherical harmonics and the B (kr) satisfy the 
n n 

differential equation 

[ ~ t ~2 + k 2 - n(n -H B (kr) = 0 . ( F-19) 

Ldr r s J n 

The general solution of equation ( F-19) can be given as a linear combin- 
ation of i (kr) and y (kr) which are the Spherical Bessel Functions of 
n n 

order n of the first and second kind respectively. 

B (kr) = C krj (kr) + D kry (kr) (F-20) 

n n J n n 'n 

where and are constants. The general solutions for f and g may 
be written as an infinite series 


l = y (c' krj _(kr) + D' kry.(kr)>J "Vtp) 


U V nm n 
n, m 


00 

g = Y (c krj (kr) + D kry (kr)) L m (e,cp). (F-21) 

L> V. nm n nm n s n 

n, m 

The constants C , D , C' , and D' may be evaluated by substituting 
nm nm nm nm 

(F-21) into (F-ll) and applying the boundary conditions (F-5). Equation 
(F-21) can be viewed as an infinite superposition of modes. 


See R. F. Harrington, Time Harmonic Electromagnetic Fields, McGraw 
Hill (1961). 

See M. Abrahamowitz and I. A. Stegun, NBS Handbook of Mathematical 
Functions, p. 437. 
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TM Modes Under Spherical Symmetry 

The TM modes are obtained by setting f 0 . It the liquid has 
spherical symmetry the boundary conditions may be satisfied by using 
only one value each for n and m in equation (F-2 1) and thus the series 
for g contains at most two non- vanishing terms: 

g (c krj (kr) + D kry (kr) L m (0 t 7,) . ( F - 

V nm n nm n s n 


Using equations (F-ll) and (F-l 6) along with (f - 22 ), the components of 
the electric and magnetic fields may be written as follows: 


■1 ; 

i x' e V : 


n(n + i) 
i x e r ? g 


io) er dr d9 


ixersinB drdcp 


l jig. 

r sin 0 Sep 


I M. 

r 30* 


F-25) 


The boundary conditions are applied by letting the container walls 
exist at r = b and the liquid-gas interface at r = a ? b (if a - 0 the 
container is full and if a = b the container is empty. ). The conditions 

pH-n and Hxn continuous at r = a and r = b imply continuity of H and 

0 

H and hence that g is continuous at r = a and r = b. This is com- 

cp 

patible with the continuity of eE»n. The condition Exn continuous implies 

that E„ and E is continuous and hence that - 7 -e is continuous at 
0 cp e 5r 6 

r = a and r = b. In summary the boundary conditions are completely 
specified by 


t 





% - 


where the coefficient of y^k^a) equation ( F-22) is zero because g must 
be finite at r = 0. (Here, = ou/€^p, o applies to the region in the gas and 
k = uu/e]I applies to the region in the liquid. ) Likewise condition (F-27) 
is equivalent to 


“ a j (k a)] = kaj (ka) + D kay (ka)l (F-30) 

e oaLono j e da L nm n nm n J 


and 


0 = ii[ C „m kb in (tt » + D „m kb ' r „ (kbl ]- 


(F*3i; 


g continuous at r = a and ~ = b 

L JL g continuous at r = a and r - b. 
e Sr ° 


( f-26) 
(F-27) 


m , 


Since the (9, rp) are independent both of radial position and fluid prop- 
erties, the condition (F-26) is equivalent to 


and 


k aj (k a) = C kaj (ka) + D kay (ka) 
o n o rm n nm n 


g(b) = c _ kb j_( kb ) + D nrn kb y n ( kb ) 

nm n nm n 


(F-28) 


(F-29) 


Equations (F-28), (F-30), and (F-31) are three independent relations in 

the eight variables, C , D , a, <u , e, u, € , and u . The inhomo- 

nm nm o o 

geneous equations (F-28) and (F-30) can be solved uniquely for C and 

nm 

D and these values are substituted in equation (F-31) which then be- 
nm 

comes a homogeneous relation in six variables a, uu , e, u> « . and u . 

o o 

We will denote this relation by 


F (u), a, «, p, e , u ) = 0 
n o o 


(F-32) 


or sometimes more simply by F^fu), etc. ) = 0. For a given set of 
values for a, e, p, c q , and (which is determined by conditions in the 

container), it can be shown that F ft plotted as a function of tu is oscillatory 
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and hence there are an infinite number of solutions to equation (F-32). 

The solution for the pth zero of equation (F-32) is called x and the 

np 

field pattern obtained by substituting x and the values for C and 
' np nm 

D into equations (F-22) and ( F - 2 5 ) is called the TM mode where 
nm mnp 


n l, 2, 3 

P = 1. 2, 3 , 

m = 0, ±1, ±2 tn, 


(The range on m comes from the properties of the spherical harmonics. ) 

Since x^ is independent of m, we see that there are a number of modes 

corresponding to the same x . This number is called the degeneracy 

np 

of x . For example, the fundamental frequency x n corresponds to 
three modes, TM 0u , TM_ Ul , and TM in and hence has degeneracy 3. 
Sometimes the first subscript is dropped and the three modes are collec- 
tively referred to as the TM n mode (which is an abuse of the term "inode"). 


We now discuss the conditions under which the resonant frequencies 

x can determine the total mass. The total mass M is a function of 
np 

three of the above variables, a, e , and e. If the resonant frequencies, 

o 

u) n p, of the modes are known, then we have the following relations in the 
five variables a, e , p o , e, and p 


0 = F^ ( u)u> etc. ) = F 1 (x xa • ®tc, ) - 
- F 2 (xj^ > etc, ) — F g (x jg, etc. ) — 


= F (x , etc. ) 
n np 


( F-33) 

, etc. ) is a relation determined by measuring 


where each of the F (x 

n np 

the resonant frequency of a TM mode. Since there are five variables 

mnp 
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it is clear that at least five different modes are necessary to completely 
determine the total mass. From the properties of the Spherical Bessel 
Functions it can be shown that each of the above relations is also inde- 
pendent; therefore; five modes are also sufficient to determine the total 
mass. If further assumptions are made, fewer modes may be sufficient. 
For example, for most liquids p as p^ ss 1 reduces the number of neces- 
sary modes to 3; if it is further assumed that e Q as 1, then the number of 
necessary modes is two; finally if in addition e is known then only one 
resonant frequency is necessary to determine the total mass. 


Alternately, it may be that the interface, r = a, is known and e 
(hence the density) is unknown; if e Q *s p Q *3 p as 1 , then the density and 
hence the total mass may be determined by a single resonant frequency. 
As a limiting case of this situation, the case a = 0 indicates a completely 
full cavity and the resonant frequencies are given by 

u 


uu 

np 


_ 2 £_ 


( F- 34) 


where u is the zero of equation (F-32) considered as a function 

np 

of the quantity kb. (The quantities u are also known as eigenvalues 

of the TM "mode". ) The measured frequency f is given by f = 

„ n P n P n P 

*r^. The calculated values for u in the case a = 0 are listed in Table 
2rr np 

1 in increasing order for the lowest ten modes. (Table 1 also includes 
results of a similar analysis for the TE modes. ) The resonant frequencies 
f^p also plotted in Table 1 are for the specific case of a 48 -cm diameter 
empty container. 


We see from Table 1 that the resonant frequencies of the lower 
order modes are widely spaced. This is primarily due to the degeneracy 
and makes it feasible to simultaneously monitor several of the lower 
order modes. 


t 
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Table F-I. 


1 

Modes 

Eigenvalues 

fc . 

™n 

u i i ■ 

2. 744 

iu, 

K 

™2 1 

U 21 = 

3. 870 

V » 

i 

TE 11 


4. 493 

£ 

f 

™31 

U 31 = 

4. 973 

£ 

te 2 i 

u' = 

21 

5. 763 

! 

™41 

U 41 = 

6. 062 


™12 

U 12 = 

6. 117 

TE 31 

U 31 = 

6. 998 

1 

™51 

U 51 = 

7. 140 

i 

™22 

U 22 = 

7. 443 


Degeneracy 


F requency 
(48 cm dia. Sphere) 

f (1 0.443 GHz 

f 21 0.766 

f' j 0. 889 

f , ' 0.984 

S, > >*> 

f,. - 1.200 

f I2 " 1 *^ 10 

f'j - 1.384 

f 51 ^ L413 

£_ = 1.472 


Example* Using Hydrogen and Nitrogen 


Equation (F-32) was solved for the four lowest order modes using 

5| « 

the FORTRAN program listed in Table 2. For given values of a, €, u. 
e Q , and the program finds the zeros of F^Or, etc.) plotted as a func- 
tion of kb where 


The p n zero is 


kb = id/JTc b. 


u = u> /pTe b. 
np np 


The computer plots the quantity »u vs. p which is essentially resonant 

np 

frequency, f , vs. total mass M. Here, 
np 


Written by A. E. Hiester, 
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Table F-II 


PROGRAM PLOT 3 

DIMENSION IFILM< 13) .ITITLEt 13) »X< 100) *Y< 100) »AL)3) *RHO( 3 ) 

DATA < IFILM-24HART HIESTER, X3A74 ) 

A(N»U.F»AL) = ( l./AL*PJ(N.AL*F*U) *YPP< N.F*U )-Al**2*PPJ ( N , AL*F*U ) *Y P ( 
1N,F*U) )/(PJ<N.F*U)*YPP(N,F*U)~PPJ(N,F*U)*YP(N,F*U) ) 
B(N.U,F»AL)=lAL**2*PJ(N f F*U)*PPJ(N»AL*F*U!-l./AL*PPJ(N,F*U)*PJ(N,A 
1L*F*U ) )/(PJ<N.F*UI*YPP(N»F*U)-PPj;N,F*Ui*YP(l\l,F»U) ) 
FUN(N»U»F»AL)=A(N*U*F»AL) # P-MN*U)-iB<N»U»F*AL)*YP(N*U) 

1 FORMAT OF10.0) 

2 FORMAT <*U*» I 2 »5H ) 

3 FORMAT ( 1H1.10X.2A8//9X.11HALPHA * UNP . 10X .*RhOBAR*/ / ) 

A FORMAT (9X*F9.5*lUX»Flu#7! 

5 FORMAT ( *0U NOT FOUND IN 10J I TERAT I ONS*// IX . 6E22 * 8 ) 

6 FORMAT ( 5HJ 1* U.12.1H ) 

P=l. 

I T I TLE I 1 ) =8H RESONAN 
I T I TLE ( 2 ) -8HT FREQUE 
I T I TLE ( 3 ) -8HNCY VS M 
I T ITLEl A) =0HASS - H2 
I T I TLE • 7 ) =8HR$9H0BAR 
I T I TLE 1 10 ) = 8H $1 

ITITLE( 11 )=8HA$9LPHA 

ITlTLE<5)-ITITLE(6)*ITITLE(8)-ITITLE<9)=ITITLEt 13)=8 h 
READ 1.<AL< I ) .1=1,3) 

READ 1 ♦ ( RHO I I ) » I * 1 1 3 ) 

CALL GRAPH(l»l»3»IFILMt0.6> 

DO 60 N= 1 * A 
I D=P 

ID-ID+10*N 

ENCODE ( 8* 2» IF ILM) ID 

ENCODE (8»6»ITITlE( 12 ) ) ID 

DO 55 1-1.3 

GO TO 17.8.9) .1 

7 LTYP-0HTP SOLID 
GO TO 95 

8 LTYP-8HTP LlO 
GO TO 95 

9 LTYP-8HNBP LIQ 

95 PRINT 3.IFILM(1) ,LTV» 

LINE-0 

DO 50 J-1.99 
RHOBAR-J 

RHOBAR-RHOBAR/100.*RHO( I ) 

F-(1.-RH0BAR/RH0( I) )**<l./3.) 

UB-7.5 
US-2. 5 

FUS-FUN(N»U8»F.AL( I ) ) 

FUB-FUNIN.UB.F.ALI I ) ) 

IT-0 

10 UM-(UB-US>/2.+US 
IT-ITfl 

IF( IT. LE. 100)00 TO 15 

PRINT 5»US»FUS»UM»FUM»UB»FUB 

STOP 
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Table F-II 


(Continued) 


15 FUM*FUN ( N *UM *F * AL ( I ) ) 

IF ( ABS< FUM) .LT..0C0U1 )G0 TO 45 

IF ( FUM.GT.O. .AND.FUS.oT .o. .OK.FUM.LT.U..AND.FUS.LT .0. )GO TO <:0 
UB*UM 
FUB«FUM 
GO TO 10 
20 US*UM 
FUS*FUM 
GO TO 10 
45 X(J)=RH0BAR 
Y t J ) = AL ( I )*UM 
L INE = L INE>1 

1 F ( L l NE »NE . 5 1 ) GO TO 5u 
L INE = 0 

PRINT 3. IFILMI 1 ) .LTYP 
50 PRINT 4 * Y I J ) »X< J) 

I F < I.NE.1IGO TO 53 

CALL LGNAPhlX,Y,99.ITITLE*!FILM) 

CALL CPGRAPH ( X t 99 ) • Y (99 ) « 1 ♦ » 1 1 +4 I 
GO TO 55 

53 CALL CLGRAPHIX.Y.99) 

CALL CPGRAPH ( X ( 99 J • Y ( 99 ) »1 » • • 1+4 > 

55 CONTINUE 

IFILM(1)«8H»9, HTP 
IF ILM ( 2 ) *8H SOLID'/ 

IF ILM ( 3 1 »8H$ 1+ TP L 

IFILM(4)-8HIQU1D'/S 

IFILMI 5)«8H1* NBP L 

I F I LM ( 6 > »8H I OU I D 

CALL COMGRAPH(.75».75»6»IFILM) 

CALL SKIPFRM 
60 CONTINUE 
STOP 
END 

FUNCTION SY(N.Z) 

Yl(Z)«-COS(Z)/Z 

Y2(Z) *-COS(Z)/Z#*2-SIN(Z)/Z 

Y3(Z)«(-3./Z**3+l./Z)*COS(Z)-3./Z**2*SINCZ> 

GO TO ( 10*20 .30t40t50>N*l 
10 SY-Y1CZ) 

RETURN 
20 SY»Y2(Z) 

RETURN 
30 SY-Y3IZ) 

RETURN 

40 SY-5./Z#Y3(Z)-Y2«Z) 

RETURN 

30 SY»7./Z»(5./Z*Y3(Z)"Y2(ZI )-Y3(Z> 

RETURN 

END 

FUNCTION SJ( N»Z I 
J1(Z»«=SIN(Z»/Z 
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Table F-II 


(Continued) 


J2tZ )=SIN<Z> /Z**2-COS(Z>/2 
J3 ( Z ) = ( 3»/Z**3-l./Z )*SIN(Z >“3./Z*»2*CUbU ) 
00 TO i 10 * 20 t 30 • 40 * 5u ) N+ 1 
10 SJ=J1(Z) 

RETURN 
20 SJ = J2 ( Z > 

RETURN 
30 SJ = J3 1 Z ) 

RETURN 

40 SJ=5./Z*J3(Z)-J2(Z) 

RETURN 

50 SJ=7./Z*!5./2*J31Z)-J2(Z) )-J3<2) 

RETURN 

END 

function pjin.zi 

FN*N 

PJ*Z*SJ(N-l»Z)-FN«SJ(NtZ ) 

RETURN 

END 

FUNCTION YPIN*Z) 

FN»N 

YP»Z*SY(N-1*Z )~FN*SY(N,Z ) 

RETURN 

END 

function ppjin»z> 
fn«n 

PPJ«(FN»(FNfl.) )/Z*SJ(N«Z) 

RETURN 

END 

FUNCTION YPPtN.Zl 
FN«N 

YPP.IFN«(FN+l.n/Z*SY<N.Z) 

RETURN 

END 
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and 


f - ^ 

np 2nb/c p. ’ 
r o o 

M = gV 

where V is the volume of the tank in cm 3 . The results may then be 
applied to spheres of any size and to any dielectric fluid. 

We have assumed that u, = u = c - 1 and plotted the results 

o o 

for three different densities corresponding to solid hydrogen, triple 
point liquid, and normal boiling point liquid; this corresponds to about 
22 percent range in density. The results for the first fcur modes are 
shown in figures Al, A2, A3, and A4. It is seen that the uncertainty in 
total mass is smaller for higher modes. Qualitatively this is because 
the field patterns are spread more uniformly throughout the cavity for 
the higher modee. For example, the uncertainty in mass vs. oru 4l (or f 41 
is less than 5 percent over most of the range. This *s to be compared 
with a density change of 22 percent indicating that the resonant mode has 
a tendency to integrate over the mass of the liquid rather than the volume. 
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APPENDIX G 


APPROXIMATE METHODS FOR AN INHOMOGENEOUS DIELECTRIC 

Introduction 

The problem is to compute resonant frequencies of a microwave 
cavity containing an inhomogeneous but isotropic dielectric. We suppose 
the cavity wall is a perfect conductor, the dielectric dissipates no power, 
and the dielectric has uniform magnetic permeability. Mathematically, we 
are dealing with the boundary value problem posed by Maxwell's equations' 
in the absence of sources and with the electric vector everywhere normal 
to the wall of the cavity. For a Fourier component Ee^ 0,4 of electric field 
with (angular) frequency <u, the boundary value problem is 

VXE = -jojuH 

7XH = joueE (1) 

E normal to boundary. 

We will use the subscripts 0, 1 to denote quantities pertaining to the 
corresponding mode for a cavity containing a uniform dielectric of permit- 
tivity e 0 or e t , respectively. We assume that the permittivity e is piece- 
wise continuous and satisfies e 0 ^ e s Also we assume that ej, - eo is 
small enough so that the set of modes for the permittivities so. e, ex are 
at most slightly different from each other in shape and can be put in one-to-one 
correspondence. For convenience we will regard e 0 as the permittivity 
of free space. 

' We use the technique, notation, and units (MKSA) of Wolfgang K. tl. 

Panofsky and Melba Phillips, Classical Electricity yind Magnetism , 
Addison-Wesley, 1955. 
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The resonant frequency uu 0 of any mode in the empty cavity ! 

may be computed from a standard solution of Maxwell's equations. * 

The resonant frequency uu of the corresponding mode when the cavity ^ 

i 

is partially filled with dielectric is lower, and we wish to estimate its 
value without further extensive calculations. We prefer estimating tech- 
niques which are insenstitive to the spatial distribution of dielectric 
material in the cavity and which do not require further computation of 

electromagnetic field strengths. ! 

{ 

j 

The first method considered is adapted from a technique of i 

approximation due to Rayleigh. Two approximations of this type are I 

formed, and it is shown that one is always at least as large as the other. , 

Then the larger is shown to be always at least as large as the true value 
of o). Finally, a refinement of the last method is described, known as 

i 

the Rayleigh- Ritz method. This first group of methods gives upper bounds 

for the true resonance frequency. The lower Rayleigh estimate has not 

been proved to be a lower bound. But the difference between the upper 

and lower Rayleigh estimates is within 10 percent of the difference between 

the empty and full cavity resonant frequencies, in the case of liquid nitrogen ; 

in a sphere in a steady uniform gravitational field (See Figure Gl). One 

should bear in mind that the upper Rayleigh estimate has been proved an 

upper bound only for the fundamental (TM on ) mode. I 

The moment methods improve on the Rayleigh methods in two ways. 

First, they provide lower bounds as well as upper bounds; and second, they 
are easier to apply to higher modes. The first order moment method yields 
a lower bound for the fundamental mode which is comparable to the lower 
Rayleigh estimate, and an upper bound which is identical to the upper Rayleigh 
estimate. Higher moments have not yet been computed, t 

4 * 

4 • 

i • 

i 




The Green's function was investigated as a tool for obtaining 
lower bounds. The resonance frequencies of the cavity are eigenvalues 
of a certain linear operator L. The Green's function is used to compute 
a norm for L -1 . Preliminary results are inconclusive as to its value, 
and careful estimates must be made with an automatic computer. 

Rayleigh Methods 
3 

The Rayleigh method of equating "potential" and "kinetic" energies 
in the perturbed field suggests the following heuristic procedure. We suppose 
that E - a E 0 for some number a independent of position and compute the 
magnetic field from eq rations (1): 


H = 


VXE 

-juHi 


aVXE 0 

-jwji 


a(- jxoUho) 

-juujj, 


(JD q 

a H 0 . 

U) 


Then we eliminate a by equating the time-average electric and magnetic 
field energies , 


/ 


2 


eE 2 



HH 2 


a 


2 




But Jz 0 Eq* = j I 


p*H 0 3 » so we conclude that uu 2 is approximately equal to 
/eo" 

( 2 ) 


(1) S = <S>Q 


/ K 


Eo s 


where K = e/eo is the dielectric constant. Notice that the ratio b of H 
to Hq is also independent of position, as a is. 


^G. Temple and W. G. Bickley, Rayleigh's Principle . Dover, 1956, pp, 1-24. 
4 

All integrals are volume integrals over the region of the cavity, unless 
otherwise specified. 
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Bat if we suppose H = bHo f° r some number ^independent of position, 
equations (1) lead to 


uoo eo 

E = b E 0 

uu e 


and the ratio of E to E 0 is now dependent on position because e is. A 
computation like that of the preceding paragraph shows that i r ' is approxi- 
mately equal to 


X = 0)0 


/e o’ 


We show (Theorem 1) that o x <• x ^ uu ^ o 0 for every mode and (Theorem 2) 
that o a < o < o for the fundamental mode. The remainder of this Appendix 

describes ways to compute lower bounds for x in the fundamental mode. 

The electric field E for any solution of (1) for a given e and x is 
an eigenfunction of the differential operator 

L = — 7 x y x 
lie 

corresponding to the eigenvalue uo s . Let S be the spherical cavity including 
the inside region and also the boundary surface. Let ts be the set of all 
vector functions on S having continuous derivatives of all orders. Then & 
is a pre-Hilbert space with respect to the inner product 

(F,G) = Jf • G , 

and we will denote by © the Hilbert space which is the (, ) - completion 

of A. It is easy to show that L maps A into 0, and so L" 1 is well-defined 

on LA. But L' 1 is bounded (Theorem 3) and LA is dense so that L' 1 has a 

unique continuous extension to © which we also denote L~\ 

5 

Now L is symmetric in the sense that 


(LF, G) = (F, LG) 
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for every F, G in A and positive in the sense that (LF, F) > 0 for every 
F in D. The same is true for Li, Lo, and the inverses of all three 
operators. Moreover we have easily 

(E 0 , LE 0 ) s (E 0 , L 0 E 0 ) 

(Ej, L' 1 Ei) s (Ei, L?Ej) . 

Now observe that 

6 1 S G ° N 

Lo - wo = — Lr*o S — ■ , 
eo v e i ^ 

so that E 0 , Ej belong to *he same eigenspace of L q . But we have agreed^ 
that E©, E, E t shall correspond to exactly the fame mode. This means that 
Ei is a scalar multiple of E 0 and 

(E j, Ei) (E 0 , E 0 ) 

(E,, L " 1 E x ) (Eo, L' x Eq) 


Theorem 1: u>i* us £ us £ tDo • 


Proof : 


(E , E ) (E , E ) 

1 1 = 1 L 

(Ei,o»{ a Ei) (Ei, Li" x Ei) 


(Ei, Ei) (Eo, Eo) (Eq, Eq) 

(Ei , L " 1 Ej) (E 0 , L ’ 1 E 0 ) (E 0 ,KL^E 0 ) 

(E 0 , E 0 ) 

■■ ■ ■ = U) ? 

(Eo, Kujo a Eo) 


*D. A. Taggart and F. W. Schott, "Ferrite-Filled Cavity Resonator, " 
Applied Science Research 25 , November 1971, page 38. 

^See the penultimate sentence of the first paragraph of the introduction. 


i 

t 





V* 


% 


-i • 


Also 


(E 0J Lq Eq) (Eq, L Eo) 


(E 0 , Eq ) 


(E 0 , E 0 ) 


(Eq, ~ Lq Eq) 

(e 0 7e 0 ) 


£ 
I >. 


f 


(E 0t ^ ujq'' E) 
(E 0 , E 0 ) 


- a 


Finally, we prove < a using the Schwartz inequality: 


(Eo, E 0 ) ;i = C r K Eo. 7 Eo )“ » " E 0 8 • '! - r 


so that 


(Eo, Eo ) 

( /K F c , /K E 0 


( 7k Eo> 7 k Eo ^ 

(Eo, Eq) 


/k 


Eo M 


Theorem 2. If a is the lowest resonant frequency, a x < a £ a 
Proof . Because a* is minimal for L and E 0 belongs to D, we have 


(E 0 , (L-a^) E 0 ) - 0 which implies a e 

Similarly, we have tjo^ maximal for L^^nd 
(El (L j- 1 - E x ) s 0 


( E o , L E 0 ) 
(E 0 . E«) 


= x* 


SO 


<u x 2 * 


(Ej, £j) 

(Ea.L ' 1 Ej) 


= a 




Remark . The result id £ uj is often called Rayleigh's principle. 
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Figures G1 through G5 show how t/(d 0 and x7x 0 vary as a 

spherical cavity is filled with liquid nitrogen in the presence of a steady 
7 

gravitational field. Preliminary data indicate that the true value lies 
between these two estimates for the fundamental mode. 

The Rayleigh- Ritz method is a refinement of the foregoing which 
also produces an upper bound for u. We will consider here only the 
case of the fundamental mode. The force of Therorem 2 is that 


0) 


2 


(F,LF) 
(F, F) 


when F is any continuous vector field on the cavity which is normal to 
the walls and twice differentiable in the interior. Theorem 2 states this 
for the case in which F is the electric field of the fundamental mode for 
an empty cavity. We now consider the case in which F is a finite sum 

N 

F= I Cn F “ 

n~l 


Where F n is the electric field of the n th mode in an empty cavity, 

II F„ || = 1, C„ is a complex number to be determined later, and N is a 

positive integer. The preceding inequality then leads to 

N N 
V’ v 


L L 


u ,2 s "-l 8=1 


C. C B tu* <F,,± F„) 
— H 


where u; B is the angular resonant frequency of the nth mode in the 
empty cavity. Then 


These results have been calculated numerically by techniques described 
in Appendix H. 
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RESONANT FREQUENCY 


TM tn MODE 



o.ioi 


<0 it SO 40 SO (0 70 10 SO too 

PERCENT FULL (NORMAL GRAVITY) 


Figure G1 

Upper and lower Rayleigh approximations to the normalized frequency 
it)/ u>o as a function of fill fraction, for ths TMqii mode using liquid 
nitrogen. 
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RESONANT FREQUENCY 



Figure G2 

Upper end lower Rayleigh approximation* to the normalized frequency 
lu/uio at a function of fill fraction, for the TMoai mode using liquid 
nitrogen. 
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RESONANT FREQUENCY 


TEw NODE 


« ii ?! ji it it u i> it 5! itt 

PERCENT FULL (NORMAL GRAVITY) 


Figure G3 

Upper and lower Rayleigh approximations to the normalized frequency 
u)/<Do as a function of fill fraction, for the TEou mode using liquid 
nitrogen. 
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RESONANT FREOUENCr 



Figure G4 

Upper end lower Reyleigh approximation* to the normalized frequency 
a)/ tco as a function of fill fractif n, for the TM 0S1 mode using liquid 
nitrogen. 
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RESONANT FREQUENCY 



Figure G5 

Upper and lower Rayleigh approximations to the normalized frequency 
u)/u>o as a function of fill fraction, for the TM 04 i mode using liquid 
nitrogen. 
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v v 1 

U)2 s min 1 C C m 2 (F — F„ ) 

' c n < K n ^ 

■=1 n= l 

where the minimum is taken over all choices of C 1( C e for which 

N 

I !c,l 2 = , 

» -1 

This is easily solved by Lagrange's method of multipliers, because the 
coefficients (F_,^ F„ ) can be computed from results already obtained for 
the lower modes in the empty cavity. We propose not only to estimate 
(u this way, but also to investigate the way this estimate depends on K. 

The Moment Method 

8 q 

The moment method was pioneered by Temple , elucidated by Kato 
and generalized to higher order by Stackgold*? It derives its name from 
the set of numbers 


m„ = 


(L" Er, , E n ) 
(E 0 , Eq) 


n= 1 , 2 , 3 , 


called moments of the operator L with respect to the vector function E 0 . 
il 

Stackgold showed that 


G. Temple, "The Theory of Rayleigh's Principle as Applied to Continuous 
Systems, "Proceedings of the Royal Society of London, vol 119 (1928), 
f . 276-293. 

9 

Tosio Kato, "On the Upper and Lower Bounds of Eigenvalues, " Journal 
of the Physical Society of Japan, vol 4 (1949), pp. 334-339. 

0 

Ivar Stackgold, The Use of Moments in Estimating the Spectrum , Depart- 
ment of the Navy, O. N. R. London Branch, August 1968, pp. 1-19, Accession 
Number A0839570. 

* Lynn Taylor Wells, Extensions, Generalization, and Clarifications of Lower 
Bound Methods Applied to Eigenvalue Problems of Continuous Elastic Systems , 
Dissertation for the PhD at Ohio State Univ. (1970), pp. 66ff. 
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JL 


■/ 


i i 

r' 




< x Sn - m n + m 2 ^ a 

b -m„ 


rru -a" 


where x is the angular resonant frequency of some empty cavity mode 
whose electric field is E 0 . (The numbers a, b will be defined shortly. 
For our purposes it will be convenient to normalize the moments as 
follows: 


M 


n - „,2n 


(L n E n , En) 


uu§ n ( E 0 , Eo! 


so that the inequality becomes 


M n - 


Man - J£ 


(4-)' 


M. 




M„ + 


Ma.-K 


M » ‘(t)’ 


The real numbers a, b are chosen as far apart as possible consistent 
with the condition that uu s be the only spectral point of L which lies strictly 
between a, b. Then for the fundamental mode the best choice of a is - * 
and the best choice of b is the angular resonant frequency of the first 
harmonic jf it were known. This leads to the simplification 


M. 




- 

n \Wo/ 


2n 


5 


for the fundamental mode, where 


* -(-fj £ ' 


Thus for the simplest estimate (n=l) the upper bound is the upper Rayleigh 
estimate 
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Consider the special case of a spherical cavity half- full of liquid 
dielectric with the fundamental mode symmetrically oriented across the 
interface. Then 

/ Eg = f Eg 

LIQUID GAS 

so that for n - 1, 2 



Then if the liquid is nitrogen (K x = 1 .44), 

= 0.8472 and M 2 = 0.7411. 

The data indicate the resonant frequencies of the lowest two modes 
have the ratio 

700 MHz 
510 MHz = 1 ‘ 
so 
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V -aZk. _i_ s LJJ . , 488 

VI ' «j 0 VM, VO. 8472 ' 1 ' 488 ' 


Hence our best possible estimate of this type is 


— = 0.847 s ( — ) 

two/ \ uuo / 


£ 0.847 - 


0.7411 - 0.7178 
(1 .488 s - 1)0. 8472 


0.847 


.0233 

1.028 


= 0.824 


to/u) 0 = 0. 920 2i (ju/uuo ^ 0 . 908 s 0.906 = uj/u> 0 . 

If Vx = 1.3 we obtain uu/uu 0 5 0 * 902 and if y 1 = 1.2, uu/ul ) 0 ^ 0. 886 (see Fig. G6) . 
Oddly enough, the upper bound does not improve as n increases from 1 to 2: 


,/Mi = .920 and a/M* =; .928. 


Calculating the lower bound corresponding to n = 2 is much more 
difficult because 


u> 8 0 ( Eq 9 E 0 ) tDo ( Eq t E 0 

involves calculating higher powers of L operating on E 0 . 

The principle disadvantage of this class of methods is that it 
requires approximate knowledge of the next eigenvalue higher than the 
one being estimated. That is, we must estimate b or Yn • But the above 
results for a spherical cavity suggest that the lower Rayleigh estimate 
is indeed a lower bound for the fundamental mode. This hypothesis should 
be investigated for a spherically symmetric distribution of dielectric 
material. 
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RESONANT FREQUENCY 



Lower bounds for (d/u^ from the first order method of moments, 
using three different values of y* . The solid lines are the 
Rayleigh upper and lower approximations. 
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Green's Function Methods 


Lemma . For any cavity V of diameter D 


/ 


d 3 R / 


D 


V 


(4 T.r R - R 


4r; 


(R-V), 


and for a spherical cavity V of diameter D, 


/ 


d 3 R ' 


(4n)^ R - R 


/ ii; 


D 

4n 


1+^2 


(Re V 


Proof . Make a change of variable 

S = R' - R. 

Since RcV, the point S = 0 always belongs to the cavity. 

Then (for any shape of cavity) as R 1 varies over the cavity, S varies 
over a region which is the cavity translated by the vector -R. So 
regardless of R the region of integration is contained in a sphere of 
radius D about -R, Hence 

D 4TT 

<; _L f s* d s C = D 
* 4n J s 3 J 4n ” 4n 


i 


d 3 s 


(4n) - sr 


V-R u u 

where W denotes solid angle. 

In the spherical case, refer to Figure G7, which is a cross- 
sectional view of the region of integration. If the R ; -origin (marked by 
the vector -R) is eccentric, from the S-origin by a distance aD, then 
the cavity will be wholly contained in a pair of hemispheres centered 
at the S-origin and having radii ( ~ -4* a ) D and bD, where a 8 + b 8 = (l/2)‘ 
Then we have 
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Figure G7 

A cross-sectional view of the region of S - integration for the Lemma. 
The R*- origin is marked by the vector -R, and is eccentric from the 
S - origin by a distance aD. 
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(^)D 


f d 3 S 1 f 
J (4Tir I s|~ £ 4n J 


£. » 1 iJ IJ " 

f A_ + _L f f 
J 4tt 4r J s“ J 4" 


i • <f *> D • i + v - ,)D ■ i 


=■£ <f +a+w 


= T^r (1 + 2 » + 2b > 

where (2a)* + (2b)* = 1. Let 2a = cos 6, 2b - sin 3 and maximize 
f (8 ) = 1 + sin 0 + cos 0. Since f / (0) = cos G - sin 0 we have the maximum 

TT 

of f (— ) = 1 + </2 and the lemma is proved. 

Theorem 3. LT 1 is bounded on LD and is extended by the operator M, 
where 


MF(R) 


ie(R / )F(R / )d a R / 


4tt R-R 


Furthermore, M is bounded and 


LT 1 ll s II Mil s 


n l/Z 


•V r /*K 3 (R / ) il 7T 1 -T T 11 d 3 R /_! /-£- 

c a L J 4n R-R | „ J c * V 4tt 


regardless of the distribution of e or the shape of the cavity. Here c 
is the speed of light and D is the greatest distance between any two points 
in the cavity, K = e/e 0 , and || K|| [ J* K** d a R / . 

Proof . First we compute LM using the identity 
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7 x7 x ■ 77 • - 7" 


and the fact that 


-7 s MI'(R) 



u>:(R / )F(R / ) 

47 

Me(R')F(R') 

4r 


v * v |r" L rT ; d " R 

[ - 4tt 6(R-r') "l d 3 R' 


- U£(R)F(R) , 

where 6 is the Dirac distribution. 

Represent tne mth Cartesian component of F by F # , set g(R,R / ) - 
4 tt| r R *Y * anC1 s ^ ow ^^f eren ^ a ^ on by subscripts, referring to primed 
coordinates. Then we have the following identity in Cartesian tensors: 


(|J€ P a g, n) » B - ( ) * b §» n + 8»n* • 


Now integrate both sides over the cavity volume in primed coordinates and 
apply the divergence theorem to the left hand side to obtain 


F,g, n d a P/ = y*(iaeF.),.g, n d 3 R / ^MC F. g, n . d 3 R' . 


The left hand side vanishes if the surface of integration is described just 
outside the (lossless) cavity. The result expressed in vector notation is 



Me F) v'g d® R' 



Me F • 


7'7'g d* R 


/ 


where 7 


/ 


is the nabla operating on primed coordinates. 


i 


i 
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Using this last result it is easy to show that 


77 • MF(R) = / U€ F ♦ 77 g 


MC F-7 ; '/ g 


■ j (v' -US F)v' g t 


and it follows that 


VxVxMF(R) = ue F 


LMF(R) - F 


:F - j (7' 1i eF|7' g 
| ^ (y/ • eF)7/ 8 • 


We see that for every G in D, 


7' *6 LG(R/ ) - v' • - v' W xG(R' ) - 0 , 

whence LMLG = LG. But the boundary value problem for L has a unique 
solution, so L has a trivial kernel and L” 1 exists. Thus 


MLG = G, 

where G is arbitrary in D. This shows that M extends L 1 . 

Since M extends L~ 1 , it will have at least as large a norm, and it 
now suffices to compute the norm of M. 


MF(R)ii, s J ue(R / )|F(R')| 


4flR-RT 11 ■ 


S /M 3 e 2 (R / ) | 


a -, 1/2 r f -, 1/2 

L d3R/ j [yi F(R / )l s d*R / 


implies the first inequality asserted for || M|l . The weaker bound on |j M|| 


comes from an estimate of 
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* >' 


4n R-ir 


(4~V J (R-R V 4^ 


which leads easily to the second bound on M ( stated in the theorem. 

Corollary . Suppose > c } for a fraction x of the \olumo of the cavity 

anc £ = for the remainder. Then the crude estimate* of the theorem 
takes the form 


Proof. 


_ ^ D" , 1 + a (K\ - 1) 

M i sT i • 

" " r- L 24 


y K? + / r- 


= [a K\ + (l-a)l^ D 3 , 

- n 4 f a K? + (1 -a) 1 
' U L 24 J ’ 


and the conclusion follows. 


Theorem 4. For the fundamental mode, 


■ n it l/2n 


for every positive integer n. 


,» l/2n 


Proof. 


2n (E. E) 


(L E, E) ’ 


but (L-E.E)s l|L’ # Eii • li E|i s || LT 0 ii • || E it * s ||m b || ||E|| 
by the Schwartz inequality and Theorem 3 above. 
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Thus 


,.2n 


U) 


L" 


M n 


and the assertion follows. 


R ema rk . Since iyn u L n h ^ is the spectral radius ; * ,1- . 

the reader may think the case n 1 sufficient. But the more general 
estimate is given in the hope that it may facilitate computation. 


Corollary . 


5 [' 


24 


,1/4 


l -Kl(Ki - 1 ) j 

Proof. Immediate from the corrollary to Theorem 


Applicatio n. For an empty spherical cavity 48 cm in diameter (any shape) 
a-0 and D = 0.48 meter, so 


x ^ \J_ZA_ x 300 x 10^ meter/sec 
2 tt 2n x 0.48 meter 


220 MHz. 


But theoretical calculations give an exact theoretical value of 54 i MHz. 
and experimental verification using liquid nitrogen at atmospheric- 
pressure gives 546 MHz. So our crude bound is quite crude. But we 
hope for a better result when 

^ 4tt |R-rH 11 

is calculated accurately as a function of R # . 

^ See Appendix F, Table I. 


1 
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A somewhat difl'rrmt ^roup of t hni qucs is associated with 
higher powers of M. 


M R F(R) 


f f -2 \ Mlii \ • • • K(R " > F ' R --- H 3 H, . . ,d a R„ 

J J 4-1 R-R. |4-| R ! - R ; | . . .4" i R~. f 


M n F (R) 


f 


K ( R i I . . K( R r . ) 1 F(R P ) 1 H 3 R. , . d 3 R B 


4" | F< - R i J , 4~ j R : - R; | ... 4” (R„ _ . , 


^ u -- 


1 2 , 


L 


f- f C ( F 'n M 3 Ri . . . d 3 R n J /2 


who re I n (K) 


/ / 


i 


* K L (Kj...K l (R 1 )d 3 R 1 ...d 3 R, 


4" j R - R i T (40* | Ri -R s | * . . . (40" | R n _ x - R n “ 
Continuing, *.ve have fo» a sphere of diameter D 


II M n F(R) ,(, • 

l 

.2 n 
c 

1 / 2 

In' (K) 

i/ d> 

If we call the last integral 

,, Fji"' , then for 




n- 1 

i,M\ s-~ 
c 

I ,/2 

TT D 3 

<K) 

1- 0 

n 

_ (U 

-J 



_1_ 


■x > ..M n lt 2n 

c 

In 4 " (K) 

rAr 

L rrD 3 



1 


/ A 

4 

4n 


r O 

1 > ° 

) l i m In 

(K) . 


n - 1 
2 


ifrU R„)d 3 R n 


1/2 


n- 1 
4n 


To estimate I n .observe that (Lemma to Theorem 3) for a spherical cavity 
of diameter D, 
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K 2 (R„ )d 3 R n ^ KID H/Z 

(4TT) 2 (R n _ 1 -R 7t f 4n V 4 / 


/• K s (R„_i )K £ (R n )d 3 R n -!d 3 R n 
7 (4tt ) 2 (R„_ i - Rn ) 2 (4n) £ (R n _ 1 -R n )^ ! 


KfD l+JZ 


d 3 R„ _ 


4n l 4 yy (4n) £ (R n _ 1 -R n r 


K 3 D , 1V2 , 

4^ V. 4 / 




4n R-Ri 


|| K(Rj )d 3 R x 


lim I n (K) £ 


KfD / 1 + v 


x a c 


k\d V i V 2 t 


C 6 \, 

(D/2V/K X V 1V2 / 


300 x 10 6 


2n “ nDyK x K 1-hJZ J 3 . 14( . 48) 1 . 2 


x 1 . 254 = 208 MHz. 


And replacing K x by 1 as an approximation for K in the estimate of I 
increases the answer by only a factor of 1.2 to the value 250 MHz. This 
suggests the crudeness of the estimate is due primarily to the way l was 
evaluated. Alternatively, set 


- /\s\ r f d 3 R!...d 3 R n 

J n (K) - J • • J (4n) 2 (R 1 -R s ) 2 . . .(4n) a (R n ^ - R n f 


and observe that for a sphere of diameter D (see Lemma). 
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n m ~ ■ 


- ^ — — 


D ; -f, / 7 \ 

I K\ “ J, (K) 


(K) 4? c 4 


1 


1 


•- D / 1V2 - 4n -1/2 4n 6 

l " J - ,KI 


n - 1 
^ 4n 


2 - 2 


-k(^) 


1/4 


lim J n 


1 

4n 


(K) 


One final approach in this vein is only heuristic at this point, 
formula for M n F(R) is an integral whose integrand is overwhelmingly 
important when the \ariables 

R , R x , R' R„ 

have all nearly the same value. This suggests the approximation 


M n F(R) a, u n eS K n (R)F(R)W n (R) 


where W. (R) =/•••/ 4 n , R-R, | 4n] r/. R S T.T. 4n | R." -R. f ' 

This strictly positive function can be estimated by the technique of the 
lemma: 


The 
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’ d 3 R 1 
4tt| R-R 1 


<r*> D 


/ |s|d|s| f -§■ +/ |s|d|s| / 


— ( - + af D 2 • - + - b 2 D s • - 
2 y 2 1 2 2 2 




7 D 2 [ 7 + a + a 3 + b s ] 
4 1 4 

i D 2 ti + * + i' 


(*) 


because a s — . 


D 2n 


Therefore W n (R) s ^ -y j 


And since 


2n (E, E) (E, E) 

“ ” (M"E, E) “ (c“^ n K n W n E,E) 


c (E, E ) 1 ^ 2n 
“ * (K“ W„ E, E) 


But for a two phase system 
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(K n W n K, E) 


1/ 


’’ r k; / 


W ri E* + 

n +s 11 | 

U^I n GAS J 


/ w.e-S-kJ / 2 


£ n 


[ K! vfy / E " 


D ' 

+ I — i 
\ 2 


2n 


L NU D 


/ 

GAS 


1 /2n 


E" 


1/2 D 

and this upper bound approac hes K x ' (^— as n increases without bound. 

Hence 


a ^ 


(D/2)7K : 


So for an empty cavity 


x 300 x 10 * 

2 tt ~ 2 tt ( 0 . 24 ) 1 . 2 


= 166M Hz . 


Whether the validity of this approach can be established is unknown, and 
even if it can be, there seems less likelihood of getting a good estimate 
for W n (R) than for J n (K). 

Recommendations for Further Study 

1. Compute the Rayleigh estimates for several of the lower modes in 
the case of a spherically symmetric distribution of liquid. 

2. Compute the bounds 


M n - 


M 2 a - M l JjL 


s M. + 


GO*-”- (;) ' "-(=>)' 


for n = 1, 2 in several of the lower modes for both normal gravity fill and 
spherically symmetric distributions of liquid. 

3. Investigate the Green's function method by computing || M|l, J„ (K), W n (K) 
accurately. 
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APPENDIX H 


NUMERICAL ANALYSIS OF THE SPHERICAL CAVITY FOR 
THE LOWEST ORDER MODES 

Electric Field Contours 

Equation (F-22) was used to generate the graphs plotted by sub- 
routine TMTEPLOT. These are plots of the electric field contours, 
|E(r, 0) | = | E | , for the TMgji, TMosi> TMoai, TMq 41 , and TEon 
modes. Choosing C Bl as imaginary (n = 1, 2, 3, 4) in (F-22) the 
equations used in TMORTE for the TM modes became 

g = z j B (z) P B (cos 0) (H- 1) 

z = u,! g- for n = 1, 2, 3, 4 

and u Bl are the eigen values of the different modes. (The outside radius, 
b, is normalized to 1. 0. P B and j„ are the Legendre and spherical 
Bessel functions, respectively. ) With 

| E | sVeFTe^Te 3 " 

1 ■ r 0 cp 

being the magnitude of the electric field, and using (H-l), the com- 
ponents become 

E r = j, (a) P, (cos 0)J 

- ^ j, (z) P, (cos 0) 


H-l 



Also 


E e = 7 777? [ 2 j * '*> p - (cos 6| ] 

* r « [ P * < COS 6) ] j z !■-> <2) - " j» (I) | 57 

and 

E = 0. 

cp 

The analygous equations for the TE modes are 
f = u,l j B (z) P n ( cos 6) 

The magnitude of the electric field becomes |e| = 

The components 

E f = E 0 = 0 

and 

E „ = u -> j- <z) ^[ p » (co * e) ] • 

The plots of | E | for the first six modes are shown in figures 8 
through 13 in the body of this report. Before plotting, |e| was scaled 
by its largest value and normalized to 10. Hence, the contours are 
numbered 1, 2 10. 

Upper and Lower Raleigh Approximations for the Normal Fill Geometry. 

The Resonant Frequency vs. Percent Full (Normal Gravity) 
graphs generated by subroutine LINPT were obtained numerically 
Integrating | E | a over different portions of the resonating spherical 
cavity for the various modes. The graphs display the upper and lower 



H-2 


,x 


l 


if 


Rayleigh approximations for the true solution. Working in spherical 
coordinates, volume is a function of (r, 0, cp ) , 

- ie, V = f (r, 0, cp). 

However, assuming symmetry in the cp direction, V = 2rr g (r, 0). 
Hence, the integrals for the different portions of the resonating cavity 
become 


where 


0. R 
k 


G E (V a ) j = 2n J f \e\ 3 r 2 si 


sin 0 drd0 


1 1 “ 

(/ 


-H/cos P 


E I* r 2 sin P dr dp 


rr 


— £ P ^ 0^ (see figure HI for 0^) 


I > 


and 


TT 

F f E ' V CP ■]■*"/ / 


E | a r 3 sin P dr dp. 


H/cos P 


G defines an integral of the field magnitude in the empty part (Vx) of the 
cavity and F defines an integral of the field magnitude in the part of the 
cavity containing fluid (V^p). Figure HI displays Vj and V^p. The 
upper and lower approximations were then computed by evaluating. 


Vi 


F + G 


and 


V3 

^ KF + 


* 

i 
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T , « 







at steps of AV__ = 1%, Here C = F + G and K is the dielectric constant 
CF 

of the fluid in the cavity. 

Figures G1-G5 show the upper and lower Raleigh approximations 
for the first five spherical modes for the normal fill geometry. 

Computer Routines 

Ten computer routines were developed and coded in Fortran to 
calculate the magnitude of the electric fields ( |E(R, 6) | ) and to integrate 
over the volume of the spherical cavity. These routines are as follows: 

TMORTE, main program, calculates |E| for the transverse 
magnetic or transverse electric modes; 

BJ, PN, and DPN are three function routines called by TMORTE 
to evaluate spherical Bessel functions of the 1st kind and 
Legendre polynomials and their derivatives, respectively; 

TMTEPLOT, subroutine called by TMORTE, plots |e| for the 
transverse magnetic and transverse electric modes; 

RCTOUR, subroutine called by TMTEPLOT, searches for 
specific contour values for TMTEPLOT to plot; 

VOLUME, subroutine called by TMORTE, inegrates |E| 3 over 
the sphere to evaluate functions of the resonating frequency 
in the empty part and the full part of the spherical cavity; 

RINTGL, subroutine called by VOLUME, performs that part of 
the integration along the radial lines, R, of the cavity; 

LINPT, subroutine called by VOLUME, plots the upper and 
lower bounds of a Rayleigh approximation to the true solution 
for a Frequency vs. Percent Full graph; 

INTRPL, subroutine called by VOLUME, interpolates between 
values of the percent full calculations to obtain required values. 
A listing of the computer routines is given as follows: 
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PROGRAM TMORTF 


DIMENSION UUI 


I FMT i ( 6 I » I FMT 2 • 6 • « VUI 


COMMON /DATA/ EMAGNI 109,181 I , I KLAB t 3 ) » ICLABI 19) .MODE, DEL TAk .DEL TATTM 
1,NRINC«NTINC.NKIP1 ,nI I P 1 .PI ,JJ( 181 ) .ESCALE TmT 

COMMON /INDEX/ INDEX. EFSWT .GAMA1 .CAMA2 .GAMA3 TMT 

tmt 

DATA (U=2. 744, 3. 870, 4. 9 73. 6. 062 > . I B* 1 • » • I P I ■ 3 « 14 159265 ) « < IRLAB*18hTMT 
1RADIAL DIRECTION I » < I CLAB ( 1 ) « 1H ) . I ICLABI 2 > ■ 1HT ) . I ICLABI 3 )■ 1HHI < I TMT 
2 ICLABI 4 • * 1HE ) ,( I CLAB <5 ) = 1HT ) , I I CLAB 1 6 ) ■ 1HA ) , I ICLABI 7 I -1H ) . ( ICLABI TMT 
38 ) * 1H ) . I ICLAB ( 9 ) = 1HD) » ( ICLAB ( 10 ) - 1HI ) . I ICLABI 1 1 ) « 1HR ) . I ICL AB I 1 2 > * TMT 
41 HE). I ICLABI 131 = 1 HC ) ,( ICLABI 14 ) = 1hT ) . I ICLABI 1 5 ) = 1H I ) . I I CL AB I 16 ) = 1HTMT 
50 ! . I ICLABI 17) = 1HN) . ( ICLABI 18 ) = 1H > * I ICLABI 191 »1H ) . I I NCR *9 ) , I I NC T » T M T 
610) * I v*4. 493. 5. 76?, 6. 998 .0.0) . I I FMT 1 =4 7H I 8H1 TM MODE »2H I.I2.1HI//5TMT 
77X.3A8/10X.13I 13.6X) ) ) . ( IFMT2*2SH( 1X.A1 .1 3.1H0.2X.13E9.2/ I ) TMT 

TMT 

INTEGER TMTt. EFSWT TMT 

REAL KX TMT 

TMT 

••••.. TMT 

ROUTINE TO calculate components OF the ELECTRIC field of the tmt 

TRANSVERSE MAGNETIC MODE UNDER SPHERICAL SYMMETRY. EQUATION TMT 

IA-22) OF NBS REPORT 9793 BECOMES TmT 

TMT 

G * P [COSI THETA)] * Z * J (2) TMT 

N N TMT 

TMT 

where z « u r/b. tmt 

N1 TMT 

TMT 

TMT 

ALSO CALCULATES COMPONENTS OF THE ELECTRIC FIELD FOR THE TMT 

TRANSVERSE ELECTRIC MODE UNDER SPHERICAL SYMMETRY. THE TmT 

ANALYGOUS EQUATION OF (A-22) THEN IS TMT 


V R/B AND B IS THE RADIUS OF THE SPHERE. 
N 1 


S EQUATION OF (A-22) THEN IS TMT 

tmt 

P [COSI THETA)] * K « J (Zl TMT 

N N N TMT 

TMT 

• K R » V R/B AND B IS THE RADIUS OF THE SPHERE. TMT 
N1 N 1 TMT 

TMT 

TmT 

tmt 

TA THE VALUE ISO DEGREES (IN RADI ANSI TMT 

TMT 

TMT 

t TMTE. MODE. NRINC.NTINC.DFLDR. EFSWT .6AMA1.GAMA2.GAMA3 TMT 

tmt 

1* CALCULATE E (AND NORMALIZE) EVERY TIME TMT 

2* CALCULATE E (AND NORMALIZE) AND WRITE E ON MAG TAPE TMT 
I LOG* UNIT NO. 1) TMT 

3. READ (NORMALIZED) E FROM MAG TAPE (LOG. UNIT NO. 1) TMT 

TMT 


R-1.0 

ASSIGN THETA The VALUE ISO DEGREES (IN RADI ANSI 
THETA-PI 

READ (60.19) TMTE .MODE . NR I NC . N T I NC . OFLDR * EFSW T .GAMA 1 • GAMA2 . GAMA 3 




t 

l f 


- I 
l f 

i 


c 

r 


9 

C 

c 


IF (EOF. 60) 18.2 
N*MODE/10 

ROT ATE « OF LOR* P 1/180.0 
IF (TMTE.EQ.2RTE) GO TO 3 
Tm MODE PARAMETERS 
INDEX.! 

L)X*U( N ) 

IFMT 1 ( 1 ) .8hl 8Hl TM M 
(30 TO A 

TE MODE PARAMETERS 
INDEX. 2 
IJX" V( N I 

I FMT 1 ( 1 I.8M8H1TE M 

KX*UX/B 

INCR-9 

INCT.JO 

IF (NRINC.NE.12 I GO To 5 
INCR-1 
I NC T ■ 1 

DELTAR.R/NR I NC 
NR1P1*NR I NC* 1 
DELTAT.THETA/NTINC 
NT IP1 *N T I NC* 1 
DO 6 1*1 .NT I Pi 
JJ( II*!-1 
CONTI NUF 

IF (EFSWT.EQ.3) GO TO 15 


i 


\ \ 
« I 


COMPONENTS OF THE ELECTRIC FIELO 

EMAX.0.0 

DO 12 I R* 1 .NR I PI 
R*( IR-1 1 *DELTAR 
Z-KX*R 
DZDR*KX 

DO 12 IT*1 .NT I P 1 

THETA. ( IT-1 I* DElTAT. ROTATE 

GO TO (7,10> » INDEX 

RADIAL (R) DIRECTION 
IF ( R.NE.0.0 ) GO TO 9 
IF CN.EO.l I GO TO 8 
EMAGNI IR . I T I *0. 0 
GO TO 12 

EMAGNI IR, IT).(2.0*KX*UX/(3.0»B) ) 

GO TO 11 

ER.-N* (N*l i* (XX**2) *PN(rOS( THETA 1 »N • *B-* ( Z »N I /Z 
ANGULAR (THETA) DIRECTION 

ETHETA»DPN(COS( THETA) ,-$IN( THE T A ) * N ) • I DZOR/R I * I Z*BJ ( Z »N- 1 ) ■ 
INI I 

EMAGNI IR. !TI»SORT(ER**2*EThETA*»2) 

GO TO 11 
TE MODE 


TmT 


Tm! 

56 

TMT 

5 1 

TMT 

5e 

TMT 

59 

TMT 

60 

TMT 

6! 

TMT 

62 

TMT 

63 

TMT 

6 k 

TMT 

65 

T m T 

66 

TMT 

67 

TMT 

68 

tmt 

69 

TmT 

^0 

TmT 

71 

tmt 

7 2 

TMT 

73 

TMT 

74 

TMT 

75 

TMT 

76 

TMT 

77 

tmt 

78 

TMT 

79 

TMT 

80 

tmt 

81 

TMT 

82 

tmt 

83 

► tmt 

8 k 

TMT 

85 

TMT 

86 

tmt 

87 

TMT 

88 

TMT 

89 

TMT 

90 

TMT 

91 

TMT 

92 

TMT 

93 

TMT 

94 

tmt 

95 

TMT 

96 

TMT 

97 

TMT 

98 

TMT 

99 

TMT 

100 

TMT 

101 

TMT 

10? 

• TMT 

103 

TMT 

104 

TMT 

105 

TMT 

106 

TMT 

107 

tmt 

108 
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K 


\ 

r 


» 



10 

c 

n 

12 

c 


n 

r 

c 


14 


c 

15 
C 

16 


17 

C 


C 

18 

C 

10 

20 

21 


C 

c 

c 

c 

c 

c 

c 


EPHI »DPN< COb< THFTA ) »-S IN< THETA ) *N» 2 *N ) 

ER*ETHETA*0«0 
MAGNI I R ♦ 1 T ) ■ ABS< EPM I ) 

IF IEMAGNI IRtin.GT.EMAX) EMAX« EMAGN U R • I T> 

continue 

PRINT IFMT1* MODEtIRLAB*(JJ(JlfJ*ltNRIPltINCR) 

00 13 I T» 1 * NT I P 1 • I NC T 

1 * I T / 10 + 1 

IF M NC T *EO* 1 ) I ■ I T 

PRINT IFMT2* ICLABI I ) * JJ( I ) • « EMAGN ( J * 1 T ) . J* 1 * NR I Pi • I NCR ) 
CONTINUE 

SCALE E ACCORDING TO LARGEST VALUE AND NORMALIZE iO 1U.C 
ESC ALF* 10*0/ EM AX 
00 14 1 R» 1 1 NR I PI 
DO 14 I T» 1 • NT I PI 

EMAGN i IR* I D «EMAGN < IR t I T > *E SCALE 

continue 

CALL TMTEPLOT 
IF IEFSWT.EU.1) GO TO 16 

RRITE <1) MEMAGNi IR»!T)*IR*l*NR!Pl>.lT-ltNTlPl) tfSCALE 
GO TO 16 

READ <1) (JtMAGNIlR#lTl.lR«l*NR|Pl)tlT«l*NTlPl).E6CALE 

PRINT IFMT 1 * MQOEt IRLABtUJl J) • J« 1 *NR JPl 1 1 NCR > 

00 17 iT.l.NTtPltlNCT 

1 ■ I T ✓ lO^l 

IF I I NCT «E0« 1 ) I« ! T 

PRINT IFMT2* I C LAB ( I ) • J J I I) • I EMAGN ( J ♦ I T ) . J« 1 *NR I P 1 . I NCR ) 

continue 

PRINT 20 

CALL VOLUME 
GO TO 1 

PRINT 21 
CALL EXIT 

FORMAT <fi2'l3»2I5»F9.0'|S«3ElO,3> 

FORMAT HMD 

FORMAT I///HM END OF JOB I 
END 

SUBROUTINE VOLUME 

SPECIAL ROUTINE TO EVALUATE I FUNCTIONS OF I VOLUMES CUT FROM 
A sphere by a plane PERPENDICULAR to The VERTICAL axis of 

SYMMETRY* HERE 

V ■ F!R*TMETA*PM!> • 2*PI *G< R • THETA »• 


TMT 109 
TMT no 

tmt ni 
tmt n? 
tm t m 

TMT in 
TMT 115 
TMT 116 
TmT 117 
TMT 118 
TMT U9 
TMT 1?0 
TMT 1?1 
TMT 122 
TMT 123 
TmT 124 
TMT 125 
TMT 126 

tmt 127 
tmt 128 

TMT 129 
TM! 130 
TMT 131 
TMT 132 
TMT 133 
TMT 136 
TMT 135 
TMT 136 
TMT 137 
TMT 138 
TMT 139 
tmt no 
tmt m 
tmt 162 
TMT 163 
TMT 166 
TMT 165 
TMT 146 
TMT 167 
TMT 168 
TMT 149 
TMT 150 
TMT 151 
TMT 152 
VOL 1 
VOL 2 
VOL 3 
VOL 6 
VOL 5 
VOL 6 
VOL 7 


VOL 8 

DIMENSION SVS$0E<1>. SV&SE'ltOl. SVICOEUl. SVICEU80I. SVcFOUIt VOL 9 
1&VCFI180I* XTRAOlt SVCF^E < 1 » • iVCFE<l80l« XTRa|< 3>» XU»1». XOOltVOL 10 







JL 


\ 


f 


I 



1 i 


1 


i 


l 

f 


i, 

I 

£ 

i 

* 


4 

$ 

I 


? 

t- 


t 



2 Y'j(i)* Y t 3 v ) t u ( 1 • * 1 l * DLOlMlt FVCFOip* fvcM 

i u ) ♦ 

r,v 1 o ( i i * Gy i i 

iV'^L 

1 1 


301* DFVCF^n* DGVKSl)* Y 1 1 1U 1 ) * Y2UU1)* YUHIIH 1 * 

) * '1LB<M lu i ) *v^L 

W 


4 Y 1 L83 1 101 ) 



VOL 

1 -3 

c 




VOL 

1 A 


COMMON /DATA/ LMAf,N( 109*1811 fiRLABn>»ICLAD( 19) 

tMODL 

• DLL TAR *Dt L T A T VOL 

1 b 


1 *NK1NC *NT I N'^ *NH ! Pi *NT I P 1 *>M * JJI 18 1 1 ♦ ESC ALE 



VOL 

16 


COMMON / I NO LX/ INDE X *EK SWT *CAMA1 * GAMA? .CAM A 3 



VOL 

1 7 

c 




VOL 

1 8 


tQoJVAlENCC < MR I NC *N I R ) . ( NR I P 1 *N 1 Hr 1 ) * INTINC* 

NIT)* 

INI I P 1 * N I T P 1 > VOL 

19 


1 * 



VOL 

2 0 


FQUlVALf MCf (VOL*TVr.cE ) 



VOL 

2 : 

e 




vci 

22 


DATA IF38*0#3?S) 



VOL 

23 

c 




VOL 

2 a 


D2R-P 1/180*0 



VOL 

2‘ 


T0X-PI 



VOL 

26 


RflX-1 .0 



VOL 

2? 


TPI38-2*0<*PI*F38 



VOL 

2? 


IPSWT-1 



V^l 

2 9 


TV0L-4.0*PI /3*0 



VOL 

3 1 

c 

RE-SCALE TO ORIGINAL EMAGN AND SQUARE 



VOL 

3 1 


rscalf-i.o/escale 



VOL 

3 : 


DO 1 IR«l*NR!Pl 



VOL 

33 


DO 1 I T * 1 .NT I P 1 



VOL 

3A 


EMAGNI 1 R * I T ) = ( EMAGN 1 I R ♦ I T ) «RSC ALE 1 **2 



VOL 

3 S 

1 

CONTINUE 



VOL 

3 6 

c 




VOL 

’ 7 

c 

SET THE INITIAL AND TERMINAL VALUES FOR P AND THETA 

VOL 

38 


RAaO'O 



VOL 

39 


RB-1.0 



VOL 

<*r. 


KtMETA«9f) 



VOL 

41 


TA»Pt /2.0 



VOL 

42 


TB«P| 



VOL 

A3 


AYZ-0.0 



VOL 

A A 

c 




VOL 

45 

c 

CAlCUlATt NECESSARY INDICES FOR A SIMPSONS 3/dTHS 

INTEGRATION 

VOL 

A6 


JB«(RA/RBXI»NIR»1.01 



VOL 

A 7 


JMAX«IRB/RBXI»NIR»0.01 



VOL 

Afl 


ib»ita/tbxi«nit*uoi 



VOL 

49 


IE»(TB/TnXl»NIT*0.01 



VOL 

50 


IEP1-IE+1 



VOL 

51 

c 




VOL 

52 

c 

SOLVE FOR THE VOLJME OF The SPHERICAL segment* U. 

LE. THETA. LE.P1VOL 

53 


I SWT-1 



VOL 

54 


DO 2 1-1*91 



VOL 

55 


CALL RINTGL II SWT * I » JB . JMAX .XY2 • V55E » 



VOL 

56 


I SWT » 3 



VOL 

A 7 


IM1-I-1 



VOL 

56 


SVSSE1 I Ml )-F38*VSSE»*INI IM1 «OEL T AT > *DEL TAT 



VOL 

5V 

2 

CC.lT 1NUE 



VOL 

60 


VALO-O.O 



VOL 

6! 


DO 3 1-1*90*3 



VOL 

62 


VALU-VALU*SVSSE< 1-1 ) *3 .0* t SVSSE t I USVSSEt 1*1 

iVSSE (1*21 

VOL 

63 

1 

continue 



VOL 

64 


I? -9 


t . 

i 


I I I 






V0L$SF»TP13B*VALlJ 


VOL 

65 


TVOLE»VrH.SS£*VOl_SS£ 


VOL 

66 


U«91 


VOL 


r 



VOL 

oft 

c 

SOLVE FOR THE REMAINING VOLUMES. PI'2.L 

.thfta.le.pi 

VOL 

b'> 


N«- 1 


VOL 

7v 


DO 10 NTMETA«KThETA«1EO 


VOL 

71 


N»N* 1 


VOL 

7 ? 


VI 3U-N > *Nl MET A 


VOL 

73 


M«-RB*C0SINTMETA*n2R ) 


VOL 

7** 


KTmP1*NThETA*1 


VOL 

7 6 


nt» :mpi 


VOL 

76 

c 



VOL 

7 7 

c 

SOLVE FOR The volumes* nThe t a.l e . THE t a. 

LL.P1 

VOL 

70 


DO 7 I-KTHPl.IEPl 


VOL 

79 


IM1-I-1 


VOL 

80 


TIM1-SINI IM1«DELTaTI»DElTAT 


VOL 

91 


IF INTHETA-90) 5. A. 5 


VOL 

ez 

4 

R«0.0 


VOL 

83 


GO TO 6 


VOL 

64 

5 

BETA* IM1*G2R 


VOL 

05 


R*-H/COS! BETA 1 


VOL 

86 

6 

JE«(R/RB*I»NIR*0.01*0.6 


VOL 

87 

C 

VOLUME OF THE SPERICAL SEGMENT CONTAING 

fuel 

VOL 

60 


CALL RlNTGL 1 A. I . JE + 1 . JMAX.VCF .VCFE 1 


VOL 

89 


T38*F38*T I M 1 


VOL 

00 


SVCFI IM1 1 *T38*VC r 


VOL 

91 


SVCFEI IM1 l«T38*VCFE 


VOL 

92 

7 

continue 


VOL 

93 

C 



VOL 

94 


VAU'2*0.0 


VOL 

95 


VALU2E-0.0 


VOL 

96 


DO 8 I>KTHP1*IF«3 


VOL 

97 


VALU2»VALU2aSVCF< |-1 » *3 .0* < SVCF ( I USVCFI 1*11 )*SVCFI 1*21 

VCL 

90 


VALU2E«VALU2E*SVCFCI 1-1)a3.0*(SVcFEI 1 IaSVCFEI 1*1 ) USVCFEI U2> 

VOL 

99 

D 

CONTINUE 


VOL 

100 


VOLCF ■ TPJ 3§*VALU2 


VOL 

101 


V0LCFE*TPI3««VALU2E 


VOL 

102 

C 



VOL 

103 


PC T F • VOCC F / T VOL 


VOL 

104 


VOU-TVOLE-VOLCFE 


VOL 

105 

c 



VOL 

106 


x:3o-ni»pctfmoo.o 


VOL 

107 


fvcfini-volcfc 


VOL 

100 


0V1 INf*V0C1 


VOL 

109 


60 :0 19.10* • JPSWT 


VOL 

no 

9 

WRITE <61*211 NThETA»VUl5»PE*V0L1.V0LCFE*PCTF.N 

VOL 

in 

10 

CONTINUE 


VOL 

112 

c 



VOL 

U3 


XU i-XU J*I.0E-0 


VOL 

114 


X<2)«Xt?>43.OE-0 


VOL 

113 


XI3I-XI3I+1O*OE~0 


VOL 

116 


DO 11 I*1'101 


VOL 

117 


UI I >•!-! 


VOL 

116 
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1 




1 




i 


i 


ii 





Y1LB1 U )=FLB1=S0RT(0L81 ) 

VOL 

1 73 

i 

17 

QLB2 = FM1- ( FM2-FM1SCJ ) ✓ < ( GAMA 2* *2- 1 • 0 ) *FM1 ) 

VOL 

174 

i 


FLB2»1 .OF+AO 

VOL 

176 

l 


IF (QLB2.LT. 0.0 ) GO TO 18 

VOL 

176 



N2*N2+ 1 

VOL 

177 

f 


Y1LB2< I )=FLB2 = SQRT(QLB2) 

VOL 

178 


18 

QLB3=FMl-( FM2-FM1SQ) / ( < GAMA3**2“1 .0 ) *FM1 > 

VOL 

179 

\ 

i 

« 

FlB3=1.0E+50 

VOL 

180 


IF (QLB3.LT. 0.0) GO TO 19 

VOL 

131 



N3 = N3 ♦ 1 

VOL 

182 



Y1LB3 ( I )*FLB3-SQRT(QLB3 ) 

VOL 

183 


19 

WRITE (61.24) FI . XDFVCF .XDGV 1 . SKFPGDC • RKFPGDC • FM2 .FLB 1 *FLB2 *F LB3 

VOL 

184 


20 

CONTINUE 

VOL 

1ft* 



WRITE (61.25) VOL 

VOL 

186 


1 C 


VOL 

187 

i 


xx«o.o 

VOL 

188 



CALL LINPT (XXfXXtXX.lUtJ.JtlOO. 0.4.0. 8*1. 0*1. XX.XX) 

VOL 

189 


f 

CALL LINPT (U*Y1*101*XX*XX.XX.XX*XX.XX.2*1.XX) 

VOL 

190 

1 

f 

CALL LINPT (U.Y2.U)l.XX»XX.XXtXX»XX*XX.2tl*XX) 

VOL 

191 



N1X*101-N 1 + 1 

VOL 

192 

t 

i 

j 

CALL LINPT (U(NlX).YlLBltNlX)tNlfXX.XX»XX*XX.XX.XX.2.2.1R.) 

VOL 

193 

i 

I 

\ 

N2X«101-N2+1 

VOL 

194 

i 

\ 

CALL LINPT ( U ( N 2 X } .Y1LB2‘N2X) »N2 . XX .XX .XX .XX • XX . XX .2 » 2 . 1R* } 

VOl 

195 

i 

i 

N3X*101-N3+1 

VOL 

196 

i 


CALL LINPT ( U ( N3X ) »Y1LB3(N3Xi .N3 . XX . XX . XX . XX • XX . XX . 2 • 2 . 1R+ ) 

VOL 

197 


} 

CALL LINPT ( XX.XX.XX.XX.XX.XX .XX. XX .XX. 3. XX. XX ) 

VOL 

19e 

j 

j k c 


VOL 

199 

i 

f 

RETURN 

VOL 

200 

j 

1 c 


VOL 

201 

j 

4 21 

FORMAT (/I4.3F10.5.2PF8.2.I5) 

VOL 

202 


i 22 

FORMAT (F8.1.F8.2) 

VOL 

203 


J 23 

FORMAT (X0H1TM MODE ( . I 2 * lH ) /3X . 3HPCT .4X .6HF ( VCF ) *5X »5HG ( Vl I . 3X ♦ 

17VOL 

204 

\ 


1HSURTMF/K ♦ G ) /C ) .3X . 19HSURT < C/ < K*F ♦ G>> • * 1 2X » 2HM2 .9X.3HLB1 • 9X « VOL 

205 

. i 

* 

23HLB2.9X.3HLB3/ t 

VOL 

206 


24 

FORMAT (F6.1.F10.3.F10.3.2F16.6.8X.F17.6.3IAX.F8.6)) 

VOL 

207 


l 25 

FORMAT ( //9H * K-1.A4/5H C*.F9.5) 

VOL 

208 


i 

FND 

VOL 

209 

* 

'1 

SUBROUTINE RINTGL ( ISWT. 1 .J9.JE.SUM.SUME1 

RIN 

1 


- c 


RIN 

2 


4 c 

ROUTINE TO INTEGRATE ALO tG THF RADIAL LINES. Rt INVOLVED 

RIN 

3 


i c 

IN A TWO DIMENSIONAL INTLSSau TO SOLVE FOR (FUNCTIONS OF) 

RIN 

4 

\ 

1 1 c 

VOLUMES OF SPHERICAL SECTORS OR SEGMENTS 

RIN 

5 

•i 

i ; c 


RIN 

6 

\ 

1 

DIMENSION RO ( 1 1 * RI108U REOU). RE 1 108 ) 

RIN 

7 

, ; c 


RIN 

8 

i 

1 

COMMON /DATA/ EMAGN 1 109 . 181 1 . I RLAfi I 3 > . ICLABl 19) .MODE *DELTAR .DELT ATR I N 

9 

j 

, 

l.NRINC.NTINC.NRIPl.NTIPi.PI.jJI 181 ) .ESCALE 

RIN 

10 

J 

i 

>! c 


RIN 

11 



GO TO (1.3.11.3). ISWT 

RIN 

12 


1 c 


RIN 

13 


i i 

DO 2 J-1.109 

RIN 

14 


i 

JM1 • J— 1 

RIN 

15 


* 

RIJMlI.I ( JM1*DELTAR|»*2)*DELTAR 

RIN 

16 


2 

CONTINUE 

RIN 

17 



I 

i 

\ 

*i 


| 

i 



s 
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xisl..: 





I 


GO TO n 
51 IM = 0 • 0 

IF (J8-JF) 4.4*11 
IJE=< JE/3)*3 
IJB=( ( JB+?)/3>*W 
I JBX=I JB 

IF njy«UBJ“l) 7,5.6 

SUM = SUM+ 1.5 # (Rt I JB ) +W ( I Jd+1 J ) 

SUM = SUM-f 1 ,5 # RII JB+1 ) +R < I JB + 2 ) 

I JBX= I JB+ 3 

IF < I JE- I JBX ) 11.11.7 

00 8 J= I JBX . I JE * 3 

SUM«SUM+R ( J-l )+3*0*(R( J)+R< J+l ) I+RCJ+?) 
CONTINUE 

IF ({JE-IJE'-1> 11*10.9 
SUM=SUM+1.5*<R( I JE + ?:+R( I JE + 1 ) i 
SUM*SUM+1.5*R< I JE+1 )+R( I JE) 


SUME fi O *0 

IF (JB-JE) 12*12.20 

I JE s ( JE /3 ) * 3 

I JB* ( I JR+2 )/3 >*3-2 

I JBX*I JR 

DO 13 J=IJB*109 

JM1 * J- 1 

REUM1 )»R( JM1 )*EMAGN( J. t I 
CONTINUE 

IF ( ( JB- I JB ) - 1 ) 16*14*15 
SUME*SUME*-1 • 5* ( RE ( I JB > + PE ( I JB+1 ) > 

SUMF SUME+1 • 5*RE ( I JB+1 >+RE( I JB+2 ) 

I JBX* I JB+ 3 

IF (IJE-UBX) 20*20*16 
DO 17 J* I JBX * I JE * 3 

SUME*SUME + RE< J-l >+3.0* (RE ( J>+RE( J + l > ) + RE ( J+2 > 

CONTINUE 

IF ( ( J£- I JE > - 1 > 20.18*19 
$UME*SUME+1*5*(RE( I JE+2>+RE( I JE + 1 > ) 

SUME*SUME+1 • 5*RE ( I JE + 1 )+RE( I JE> 

RETURN 

END 

SUBROUTINE linpt(xi *Y I *NPT*iNXD*FMNX*FMXX.NYD*FMW*FMXY*N*LS 


LINEAR PLOTTING ROUTINE 

COMMON /DD/ IN*I0R*IT*IS*IC*ICC*IX*IY 
COMMON /DOC/ LU.LUC.IFL 

COMMON /LABEL/ LABELX I 5 ) *LABELY( 5 > *LSWTX*LSWTY 
COMMON /DATA/ EMAGN( 109*181 )* IRLABI 3 »* ICLAB( 19) »MODE*DELTAR 
1 tNRINC *NT INCtNRlPl .NTIP1 *PI *JJ( 101 I *ESCAL£ 

COMMON /INDEX/ INDEX *1 FSWT *GAMA1 *GAMA2 *GAMA3 


r i n i a 

RIN 19 
R I N 7 0 
WIN 21 
R I N 22 
RIN 2 3 
R I N 2 4 
R'IN 2 5 
R I N 26 
R I N 2 7 
RIN 28 
R I N 29 


R I N 3 2 
RIN 33 
RIN 34 
RIN 35 
RIN 36 
RIN 37 
RIN 38 
RIN 39 
'RIN 40 
RIN 41 
RIN 42 
RIN 43 
RIN 44 
RIN 45 
RIN 46 
RIN 47 
RIN 48 
RIN 49 
RIN 50 
RIN 51 
RIN 52 
RIN 53 
RIN 54 
RIN 55 
RIN 56 
RIN 57 
RIN 58 
RIN 59 
RIN 60 

* I SYCH ) L I N 1 

LIN 2 
LIN 3 
l;n 4 
L 1 1 5 

LIN 6 
LIN 7 
LIN 8 

• DELTAUIN 9 

LIN 10 
LIN 11 
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DIMENSION XIINPT)* Y I ( NPT ) * ID(4) 


DATA iID*32H RG PETERSON X3184 ) 

DATA t IXMINMOO) •! IYMIN»100> t(LSWTX*2> . (L$WTY = 2> 

DATA (LABELX = 5(8H ) ) t ( L ABF.t Y * 5 ( 8H )) 

DATA (LABELX*40N PERCENT FULL ( NORMAL GRAVITY! 

DATA ( LABEL Y*40H RESONANT FREQUENCY 

integer efs*t 

GO TO (1*10*21) * N 

IOR-O 
I N*0 
IC*0 

IF (EFSWl •EG#3) CALL DDINIT <4tID> 

I T*0 

DRAW A BOX AROUND PLOTTING AREA 
CALL DDBOX 10*1023*0*1023) 

SPECIAL LABEL PLOTTED 
I T* 1 
15*3 
I X*675 
IY*800 

LABEL-8HTM MODE 

IF ( INDEX. EG. 2> LABEL-8HTE MODE 
CALL DDTAB 

CALL DDTABNA8 ( 1 *LAB£L • 1 ) 

MODEXMOOO+MODE 

ENCODE ( 8*23 *LBMODE ) MODEX 

IS- 1 

IX-71 5 

IY-790 

CALL DDTAB 

call DDTABNA8 (1*LBM0DE*1) 

NY*NYD+1 

LTYDIV-900/NYD 

iymax«iymin+ltydiv*nyd 

NX-NXD+1 

LTXDIV-900/NXD 

ixmax*ixmin+ltxdiv*nxd 

DRAW DIVISIONS OR TICKMARKS FOR X-AXIS 

iy-iymax 

IX-MINX-IXMIN 
CALL DDBP 
IY-IYMIN 
CALL DDVC 
DO 4 I *2 »NX 
1 Y* I YMAX 

ixs*ix*ixmin^u-1)*ltxdiv 

CALL DDBP 


LIN 

12 

L I N 

n 

LIN 

1 A 

L I N 

is 

L I N 

16 

L I N 

1 7 

LIN 

16 

LIN 

19 

LIN 

20 

L I N 

2 1 

L I N 

22 

L I N 

23 

LIN 

24 

LIN 

25 

LIN 

26 

LIN 

27 

LIN 

28 

LIN 

29 

LIN 

30 

LIN 

31 

LIN 

32 

LIN 

33 

LIN 

34 

LIN 

35 

LIN 

36 

LIN 

37 

LIN 

38 

LIN 

39 

LIN 

40 

UN 

41 

LIN 

42 

LIN 

43 

LIN 

44 

LIN 

4 5 

LIN 

46 

LIN 

47 

LIN 

46 

LIN 

49 

LIN 

50 

LIN 

51 

LIN 

52 

LIN 

53 

LIN 

54 

LIN 

55 

LIN 

56 

LIN 

57 

LIN 

58 

LIN 

59 

LIN 

60 

LIN 

61 

LIN 

62 

LIN 

63 

LIN 

64 

LIN 

65 
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I Y5= I Y = I YMAX-8 

call ddvc 
ixsjxs 
I Y* IYS 

IF ( I .EQ.NX > GO TO A 
GC TO (3*2) t LSWTX 
? 1 Y= ! YM l n+8 

call DDBP 
3 I Y = I Y M I N 

CALL DDVC 

A continue 

c 

C DRAW DIVISIONS OR TICKMaKKS FOR Y-AXIS 

I X * I XMAX 
IY=MINY=IYMIN 

call ddb° 

IX= IXMIN 
CALL DDVC 
DO 7 1 = 2 *NY 
I X= I XMAX 

IYS=IY = IYMlN+t I-1)*LTYDIV 

CALL DDBP 

1 XS= I X= I XMAX-8 

CALL DDVC 

IX=IXS 

IY=IYS 

IF ( UEU.NY) GO TO 6 
GO TO <6*5> * LSWTY 

5 I X= I XM I N+8 
CALL DDBP 

6 I X= I XM IN 
CALL DDVC 

7 CONTINUE 
C 

C NUMBER THE X-AXIS 

IS=1 
I T * 0 
I0R*0 

F INCX* ( FMXX-FMNX I /NXD 
DO 8 I * 1 tNX 
IYS'IY«IYMIN-15 
I Ml ■ I - 1 

I XS* I X* I XMAX- IM1*LTXDI V-32 
FNUMB*FMXX-IM1*FINCX 
ENCODE (0.2^.FNUMBX) FNUMB 
CALL DDTAB 

CALL DDTABNA8 ( 1 tFNUMBX * 1 ) 

8 CONTINUE 

C NUMBER THE Y AXIS 

IS*1 
IOR-0 

FINCY« ( FMXY-FMNY) /NYD 
DO 9 I * 1 *NY 
I X«60 


LIN 66 
LlN 67 
LIN 66 
LIN 6 9 
LlN 7') 
LIN -> 1 
l I N 7? 
L I N 73 
LIN 74 
l I N 76 
L I N 76 
LlN 77 
UN 7e 
LIN 79 
LIN 8 J 
LIN 81 
LIN a; 
LIN 83 
LIN 8 A 
LIN 86 
LIN 86 
LIN 87 
LIN 88 
LIN 89 
LIN 90 
LIN 91 
LIN 92 
LIN 93 
LIN 94 
LIN 95 
LIN 96 
LIN 97 
LIN 98 
LIN 99 
LIN 100 
LIN 101 
LIN 102 
LIN 103 
LIN 104 
LIN 105 
LIN 106 
LIN 107 
LIN 108 
LIN 109 
LIN 110 
LIN HI 
LIN 112 
LIN 113 
LIN 114 
LIN 115 
LIN 116 
LIN 117 
LIN 118 
LIN 119 


i 

8 

\ 
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LIN 

U > 


I Y = I YMAX- IM1 # LTYDIV 

LIN 

12 1 


FNUMB=FMXY- JM1*FINCY 

L IN 

122 


ENCODE ( B *25 ♦FNUMBV ) F.NUMB 

L I N 

12’ 


call ddtar 

L IN 

] 2 a 


CALL DDTAriNAa (ltTNUMBYtl) 

LIN 

1 2 5 

Q 

CONTI NUE 

LIN 

126 

r 


LIN 

1 2 7 

r 

LABEL X AND Y AXES 

LIN 

128 


I T= X 

LIN 

1 2 9 


I S= 3 

LIN 

n j 


I X = 40 

L I N 

131 


I Y= 30 

LIN 

132 


CALL ddtab 

L I N 

133 


call DDTABNA8 ( 5 *L ABELX { 1 ) * 1 ) 

LIN 

1 3 A 


X 

ii 

o 

LIN 

135 


I Y = 40 

LIN 

136 


IOR = 1 

LIN 

137 


call DDTAB 

LIN 

138 


CALL DDTABNA0 { 5 iLABEL Y ( 1 > * 1 ) 

LIN 

139 

c 


LIN 

1 A 0 


XM I N=FMNX 

LIN 

1A1 


XMAX = FMXX 

LIN 

1A2 


YM I N=FMN Y 

LIN 

1 A 3 


YMAX^FMXY 

LIN 

1 A A 


FM I NX = I XM I N 

LIN 

1 A 5 


FmAXX= ixmax 

L IN 

1A6 


FMINY*IYMIN 

LIN 

1 A 7 


FMAX Y“ I YMAX 

LIN 

1A8 


IS = 0 

LIN 

1A9 


I T ■ 0 

LIN 

150 


I OR = 0 

LIN 

151 


GO TO 22 

LIN 

152 

c 


LIN 

153 

10 

ixl*iyl=o 

LIN 

15A 


JSWT= 1 

LIN 

155 


GO TO ( 11 *12 ) t LS 

LIN 

156 

c 

PLOT VECTORS 

LIN 

157 

11 

CALL DDCONVEC 

LiN 

158 


GO TO 13 

LIN 

159 

c 

PLOT SYMBOLS OR CHARACTERS 

LIN 

160 

1? 

I CC* I 5 YCH 

LIN 

161 


CALL DDSYMBOL 

LIN 

162 

n 

DO 20 I«1*NPT 

LIN 

163 


X*XI ( I ) 

LIN 

16A 


y«yi ( n 

LIN 

165 


I X* < X-XM I N ) / ( XMAX-XM IN ) * ( FMAXX-FM I NX I ♦FM I NX 

LIN 

166 


IY* « Y-YMIN)/ UMAX-YMINHMFMAXY-FMINY1+FMINY 

LIN 

167 


SLOPE* ( IY-I YL > ✓ C IX-IXL) 

LIN 

168 


IF UX.GE.IXMIN) GO TO 14 

LIN 

169 


I Y*SLOPE # ( I XMI N- 1 XL ) +1 YL 

LIN 

170 


IX-IXMIN 

LIN 

171 


JSWT *2 

LIN 

172 


GO TO 15 

LIN 

173 
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IF ( I X * LF • I XMAX ) G 0 TO 1 6 
I ys SLOPE# ( I XMAY- I XL U I Yl 
I X = I XMAX 
JSW T = 2 

IF ( I Y # OF • I YM I N ) GO TO 16 

I Y = I YM I N 
GO TO 17 

IF ( I Y • LF • I YMAX ) GO TO 18 
I Y = I YM A X 

i x = t iy-iylj/slope+ixl 
JSWT = 2 

IF { l.EQ. l.AND. JS wT.lU. 2 I 00 TO 19 
IF ( IX.E0.1XL.AN0. I Y.tQ. IY l > 00 TO 
CALL DDXY 
I XL = IX 
I Y|_ * I Y 
CONTINUE 

CALL DDTAB 
GO TO 22 

FRAME ADVANCE 
CALL DDFR 

RETURN 

FORMAT (13) 

FORMAT 1F4.0) 

FORMAT (F4.2) 

END 

SUBROUTINE TMTEPLOT 


LIN 17a 
LIN 175 
L I N 176 
LIN 177 
LIN 173 
LIN 179 
LIN 180 
LIN 181 
LIN 182 
LIN 183 
LIN 104 
LIN 185 
LIN 186 
LIN 187 
LIN 188 
LIN 189 
LIN 190 
LIN 191 
LIN 192 
L IN 193 
l N 194 
LIN 195 
LIN 196 
LIN 197 
LIN 198 
LIN 199 
LIN 200 
LIN 201 
LIN 202 
LIN 203- 


SPECIAL PLOTTING ROUTINE TO PLOT MAGNITUDES OF ELECTRIC FIELDS TMP 

TMP 

COMMON /DD/ IN*IORtlT*lStIC*ICCt!X.lY TMP 

COMMON /DDC/ LUtLUCtlFL TMP 

COMMON /DATA/ EMAGN( 109.101 ) • IRLAB<3> • ICLAB< 19) .MODE .DELTAR tDELTATTMP 
ltNRINC.NTINCfNRlPl.NTlPl tPI *JJ( 181 ) •CSCALE TMP 

COMMON /INDEX/ INDEX *EFSWT *G-»MA1 » UAMA2 *GAMA3 TMP 

COMMON /SUB/ IRRI 181 ) • IXS( 181 ) • I YSU81 ) *PID2*PIDNT I tR.NP TMP 

TMP 

DIMENSION LABEL15) * IXNUO). RAD(36>. 1D<4> TMP 

TMP 

DATA (ID«32HRG PETERSON X3184 ).<LABEL«40H R F MASS TMP 

1GUAGING-PL0T OF E t TM MODE ).( RAD«20(475.0 ) .477.0.478.0 .48Q.0.4TMP 
281 . 0 *483 .0 *<*85 .0.468* 0. 490. Ot 488. 0*485 .0. 483.0 .481 .0.480.0 .478 .Q.4TMP 
377.0*475.0) TMP 


IOR-0 
IC-0 
I T »0 

CALL ddinit 


(4* ID) 


BOX AROUND PLOTTING AREA 
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1 N = 1 

CALL DDBOX 10. 1023. 0.1023) 

RLuT A CONTINUOUS CIRCLE CtNTEKED AT 
WITH RADIUS = 450 RASTER UNITS 
IS = 0 

CALL DDCONVEC 
PIDNT I =PI /NT INC 
PI02=PI/2.0 
DC 1 1*1 .NT 1 P 1 
T = ( I - 1 >*P[DNT I 
TT=P I D2-T 

IX=C0S(TTI»450. 0+512.0 

IY=SIN<TT)*450. 0+532.0 

IXSINT IP1-I + 1 ) * 1 024- I X 

IYSINT IPl-I + 1 ) * I Y 

CALL DOXY 

CONTINUE 

DO 2 I *2 .NT I P 1 

IX-IXSI II 

IY*IYS< I) 

CALL DOXY 
CONTINUE 

NUMBER CIRCLE EVERY 10 DEGREES 
DO 3 1*1,36 
I DEG* ( 1-1 1*10 
TT*P I D2- I DEG*P I DNT I 
COSTT«COS(TT) 

SINTT«SIN<TT> 

CALL DDCONVEC 
IX»COSTT*450. 0+512.0 
IY*SINTT*450. 0+532.0 
CALL DOXY 

IX*C0STT*458. 0+512.0 
IY*SINTT*458. 0+532.0 
CALL DOXY 

1X«C0STT*RAD( I 1+512.0 
I Y*S I NTT*RAD ( 11+532.0 
IDEGX*1000+IDEG 
ENCODE (8.13.NVI IDEGX 
CALL DDTAB 

CALL DDTABNAS (l.NV.l) 

CONTINUE 

LABEL PLOTTED 
I T* 1 
ISO 
IX*10 
I Y»25 

LABEL I 5 • *8hM MODE 

IF ( INDEX. EO. 21 LABEL! 5 I *8HE MODE 

CALL DDTAB 

CALL DDTABNA8 ( 5 »LABEL < 1 » * 1 > 

MODEX* lOOO+MODE 


TMP 2 5 
TMP 26 
TMP 2 7 

RASTER COOK. (512.562) TMP 28 

TMP 24 
TMP 30 
TMP 31 
TMP 32 
TMP 33 
TMP 34 
T mP 36 

TMP 36 
TMP 37 
TMP 38 
TMP 39 
TMP 40 
TMP 41 
TMP 42 
TMP 43 
TMP 44 
TMP 45 
TMP 46 
TMP 47 
TMP 48 
TMP 49 
TMP 50 
TMP 51 
TMP 52 
TMP 53 
TMP 54 
TMP 55 
TMP 56 
TMP 57 
TMP 58 
TMP 59 
TMP 60 
TMP 61 
TMP 62 
TMP 63 
TMP 6. 
TMP o5 
TMP 66 
TMP 67 
TMP 68 
TMP 69 
TMP 70 
TMP 71 
TMP 72 
TMP 73 
TMP 74 
TMP 75 
TMP 76 
TMP 77 
TMP 78 
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FNCODE < 8 . 13 .LtfMODE ) MODE X 
I S* 1 
I X*79 4 
IV* 15 

call ddtab 

call DDTABNA8 ( 1 * LBMODE • 1 ) 

C 

I S * 0 

NTlD2=NTlNC/2+l 

r 

C SF ARC H FOR CONTOURS 

DO 12 NC = 1 * 1 0 
R = NC 

C ASSIGN BEGINNING INjICFS. 2 FOR 

DO 5 ITH*1*NTIP1 
IF (M0DE.EQ.41 ) GO TC 4 
IRR< I T H I * 2 
GO TO 5 

4 I PR ( I T H ) = 2 0 

5 CONTINUE 
C 

ITH«1 

6 CALL rctour 
ithb* I Th 

c 

c plot the contours 

IF ( NP *FQ*0 * GO TO 10 
CALL ddconvec 
DO 7 I s 1 »NP 

ix*ixs< n 

I Y» I YS ( I ) 

CALL DDXY 

7 CONTINUE 
IXR«1024-IXS(NP) 

IF ( IXSil J.EO.IXRJ GO TO 8 
CALL DD8P 
CALL DDCONVEC 

8 DO 9 1*1. NP 
NPMIP1*NP-1+1 
IX-1024-IXSINPMIP1 ) 

!Y«IYS<NPMIP1 I 
CALL DDXY 

9 CONTINUE 
CALL DDBP 

C 

10 DO 11 I TH* I THB tNT I P 1 

IF (IRR(ITH).LE.NRIPl) GO TO 6 

11 CONTINUE 

12 CONTINUE 
C 

c frame advance and return 

CALL DDFR 

c 

RETURN 


T M P 7 9 
TVP 8. 
T **P 8 1 

T"P 82 
T vP 81 
T^P 8 4 
T"P 85 
TVP 86 
TMP 8 7 
T^P 88 
TMP 89 
TMP QC 
T v p q i 

TM OR TE 11.21.31. AND 2U FOR 41 TMP 92 

TMP g i 
TMP 0 4 
T Mp g*> 
TMP 96 
T^P 9 7 
TMP 9 0 
TMP 99 
TMP 2 
TMP 101 
TMP 10? 
TMP 10-' 
TMP 104 
TMP 10S 
TMP 106 
TMP 107 
TMP 108 
TMP 109 
TMP 110 
TMP HI 
TMP 112 
» TMP 113 
TMP 114 
TMP 115 
TMP 116 
TMP 117 
TMP 118 
TMP 119 
IMP 120 
TMP 121 
TMP 122 
TMP 123 
TMP 124 
TMP 125 
TMP 126 
TMP 127 
TMP 128 
TMP 129 
TMP 130 
TMP 131 
.MP 132 
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13 


TMP 133 

FORMAT (|31 TMP 134 

END TMP 1 3 *>- 

SUBROUTINE RC T OUR RC T 1 

RCT 2 

SPECIAL ROUTINE To SEARCH FOR ALL VALUFS BELONGING TO A CONTOURRCT 3 

RCT 4 

COMMON /data/ tMAGM lus.ian . IRLAPI 3 1 1 ICLABI IV) fMOOt .OLLTAK.OELTATRCT 5 
l.NHINC tNTINC.NKIPl.NT 1^1 ,PI .JJI 181 1 .ESCALE KC T 6 

COMMON /SOB/ IRN< 181 1 1 I XS< 181 I * 1 YS< 181 I •PID«'«P IONT I *R*NP RC T T 


NP-0 

00 A I TH« 1 «NT I Pi 
IRB*IRR< I Th ) 

IF t IRB.GT.NRIPl ) GO TO 4 
EMI «£MA6N ( IRB-1 tITM) 

DO 1 I R* I RB *NR I P 1 
E*EMAGN ( IR.lTHl 

IF (EM1.LE.R.AN0.R.LE.E) GO TO 2 

IF (EM1.GE.R.AND.R.GE.E i GO TO 2 

EM1-E 

CONTINUE 

RNUM-0.0 

GO TO 3 

RNUM*EM1 -R 

RDEN«EM1-E 

R|NTP« ( IR-2.0*RNUM/R0EN)*DELTAR 

RADIUS«RINTP*450.0 

TT-PID2-C ITH-1 )»PIDNTI 

NP»NP*1 

IKS(NP)«C0S(TTI»RADIuS+S12«0 
IVS(NP)«SINITT)»RADIUS*-«>32.0 
IRRI ITHl.IR+1 

IF (RNUM.EQ.O.O.AND.NP.NE.O 1 GO TO 5 
CONT I NUE 

RETURN 

END 

SUBROUTINE INTRPLt I U*L *X * Y «N »U* V ) 

INTERPOLATION OF A SINGLE-VALUED FUNCTION 

THIS SUBROUTINE INTERPOLATES* FROM VALJES OF THE FUNCTION 
GIVEN AS ORDINATES OF INPUT DATA POINTS IN AN X-Y PLANE 
AND FOR A GIVEN SET OF X VALUES (ABSCISSAS I* THE VALUES OF 
A SINGLE-VALUED FUNCTION Y « Y(X). 


The INPUT PARAMETERS ARE 

IU ■ LOGICAL UNIT NUMBER OF STANDARD OUTPUT UNIT 
L » NUMBER OF INPUT DATA POINTS 
I MUST BE 2 OR GREATER) 

X • ARRAY OF DIMENSION L STORING THE X VALUES 
(ABSCISSAS) OF INPUT DATA POINTS 


RCT 8 
RCT V 
RCT 10 
RCT 1 l 
RCT 12 
RCT 13 
RCT 14 
RCT 15 
RCT 16 
RCT 17 
RCT 18 
RCT 19 
RCT 20 
RCT 21 
RCT 22 
RCT 23 
RCT 24 
RCT 25 
RCT 26 
RCT 27 
RCT 21 
RCT 2« 
RCT 30 
RCT 31 
RCT 32 
RCT 33 
RCT 34 
RCT 35 - 

0001 

0002 

0003 

0004 

0005 

0006 

0007 

OOOB 

0009 

0010 

0011 

0012 

0013 

0014 

0015 

0016 


I 
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( IN ASCENDING ORDER ) 

ARRAY uF DIMfNSIwN L STORING THf Y VALUES 

(ORDINATES) OF I NPl I T DATA POINTS 

NUMRFK UE PUINTS AT WHICH INTERPOLATION jf 

Y VALUt. (ORDINATE) I b DESIRED 

(MUST bE 1 OP GREATER) 

ARRAY OF DIMENSION iN STORING THE X VALUE 3 
(ABSCISSAS) OF OF SIRED POINTS 


Thf 


" ThF OUTPUT PARAMETER JS 


ARRAY OF DIMENSION N WHERE THE INTERPOLATED Y 
VALUES (ORDlNATcb) ARl To bE DISPLAYED 


declaration statements 

DIMENSION X ( 1 I * Y ( 1 ) .till > * V < 1 ) 

EQUIVALENCE ( PO *X3 ) * ( G* ♦ Y 3 ) • IU1 *T 3 ) » ( UK. ♦ DX ) * i I MN * X2 * A 1 « FM 1 ) * 
1 I IMXtXS * A5*FMS ) • ( J» jfc «SA U ( Y2 *W2 *W4*Q2 I • ( Y5 t *3 *Q3 ) 

PRELIMINARY PROCESSING 

10 L0*L 
LM1-L0-1 
LM2 *LM 1-1 
LP1*L0+1 
N0*N 

IF ( LM? ) 00 1 100 #100 
ion IFINO>91*91*105 
ins DO 11 !«2*L0 

I F ( X < 1-1 >-X< I ) ) 11 *95*96 

11 CONTINUE 
IPV-0 

MAIN DO-LOOP 

DO 80 K * 1 *N0 
UK-UlK) 


ROUTINE TO LOCATE THE DESIRED POINT 


20 

IFILM2>200.27.200 

200 

ifiuk - 

X 1 LO * 1205 *26 t?6 

205 

IF (UK - 

XII) >2$t206«208 

208 

IMN«2 



imx.lo 


2 1 

1 ■ 1 !MN*lMX>/2 


ifiuk - 

XIII >22»23»2T 

22 

!MX«I 



GO TO 2fc 


23 

I MN- 1*1 


24 

1 F 1 | MX - 

1 MN ) 240 » 240 1 2 1 

240 

1 * IMX 



w *17 

or 18 

y i j 

1 1 r / ( 

r t ‘ 2 3 

Im j?t* 

* ; ’<?5 
un2 6 

«)( 2 7 
00 2 B 
0 J 2 u 

uU3o 
' *31 
00 32 
) 33 
(JO 1a 
0 J 3 S 

0036 

0037 
U0 38 
U j 39 
U J A U 

A 1 
' 1,1 A 2 
J.,43 
0 U A A 
00 AS 
0 0 A 6 
00 A 7 
00 A0 
00A9 

0050 

0051 

0052 

0053 
005A 

0055 

0056 

0057 

0058 

0059 

0060 
0061 
0062 

0063 

0064 

0065 

0066 

0067 

0068 

0069 

0070 
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GO TO 30 
2b 1*1 

GO TO 30 
26 I=LP1 
GO TO 3n 
?7 1 = 2 

r 

r check ir ! = ipv 

30 I F ( [ - lPV)300t70#300 

300 IPV* I 

ROUTINES TO PIC* UP NECESSARY X AND Y VALUES AND 
TO ESTIMATE THEM IF NECESSARY 


40 J* I 

I F ( J-l J 40 1 *400 * 40 1 

400 J*2 

401 IF( J-LPl)403t402*403 

402 J*L0 

403 X 3* X (J-l ) 

Y3»YIJ-n 

X4»X ( J> 

Y4*Y ij) 

A3«X4-X3 

FM3«I Y4-Y3)/A3 

I F I L M2 ) 404 • 4 3 *404 

404 |F<J - 2 1 40 5*41*405 

405 X2-XJJ-2) 

Y2«Y< J-2 ) 

A2«X3-X2 
FM2*(Y3-Y2)/A2 
IF(J - LO ) 4 1 *42 *4 1 

41 X5*XIJ41) 

Y5*Yt J*l) 

A4»X5-X4 

FM4 » < Y5-Y4 ) /A4 
I F I J - 2 * 45 *4 10 #45 

410 FM2 *Fm3+FM3-FM4 
GO TO 45 

42 FM4«FM3*FM3-FM2 
GO TO 45 

43 FM2«FM3 
FM4*FM3 

45 IF I J - 3 1 46*46*450 
450 Al«X2-XJJ-3> 

FM1.IV2-YI J-3> »/Al 
GO TO 47 

46 FM1»FM*4FM2-FM3 

47 IFIJ - LM1 1470*48.48 

470 A5«XU^2)-X5 

FM5*I Y< J^2)-Y3)/A5 
GO TO 50 



( >r 7 1 
0/ 72 
Go 7 4 

,<> /i* 
ii f *■> 
u7fc 
r ? 7 
I , '8 
L- >11 

00 3 0 

0U81 
C~82 
0083 
u 0 8 4 

0 0 8b 

V )8fe 

OH87 

joee 

0^89 

jOQ'j 

0 ( J 9 1 
0092 
OC93 

O', ^4 

0095 

0096 

0097 
Q09e 
00 99 
0100 
0101 
0102 
0103 
0 104 

0105 
0 106 
0107 

0106 

0109 

0110 
0111 
0112 

0113 

0114 

0115 

0116 

0117 

0118 

0119 

0120 
0121 
0 122 

0123 

0124 
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NUMERICAL DIFFERENT [ATION 

0 126 
1 2 7 

AO IF ( I-LP1 ) 600 * 52 * 500 

'128 

500 W2*AHSI FM4-FM3 > 

^12R 

W 3 ■ ARS ( F M 2 " F M 1 ) 

Jl 3u 

5w»w2<**3 

v 1 3 1 

IF(Sw»61*601*5l 

v ‘ 32 

6nl w2*o#6 

- 1; ) 

W 3 * 0 • A 

* 1 ^4 

SW-1 .0 

0 1 3 

A 1 T3*lW2 # FM2+W3*FM3 ) /$W 

0 136 

IFi l-l 152*54 t 02 

013 7 

62 W3*AB$I FM5-FM4 ) 

o 1 3ft 

W4-ABSIFM3-FM? ) 

c 139 

5* a w3+W4 

0 140 

IF (54)53*520*53 

r* i l 

520 43*0.6 

0142 

W4»0.‘ 

0U3 

SW* 1*0 

u l uu 

53 T4* ( W3*FM3+*4#FM4 ) /SW 

014*5 

I F ( I -LPl ) 60 • 5 30 .60 

0 146 

*30 T3*T4 

0 147 

SA* A2 >A3 

oue 

T4*0.5*(FM4*FM5-A2* ' A2“A3)*< FM2-FM3 >/(SA*$A) > 

0 149 

X3*X4 

0150 

V3*Y4 

CM5 1 

A3*A2 

0152 

FM3«FM4 

0153 

GO TO 60 

0154 

54 T4*T 3 

0155 

SA-A3+A4 

0156 

T3«0.5«IFM14FM2-A4#IA3-A4HMFM3-FM4)/(SA*5A) ) 

0157 

X3-X3-A4 

0158 

Y3»Y3-FM2«A4 

0159 

A3»A4 

016 0 

FM3-FM2 

0161 

0162 

DETERMINATION of the COEFFICIENTS 

0163 

0164 

60 Q2- I2.0»IFM3-T3) ^FM3-T4)/A3 

0165 

Q3«I-FM3-FM3*T‘UT4 )/( A3*A3 ) 

0166 

0167 

COMPUTATION OF THE POLYNOMIAL 

0168 

0169 

70 t)X*UK-P0 

0170 

PO V IK )»0040X*( 01 4DX*( Q2^0X-Q3) ) 

0171 

RETURN 

0172 

0173 

ERROR EXIT 

0174 

0175 

90 WRITE ( IU *2090 I 

0176 

GO TO 99 

0177 

91 WRITE < IU • 209 1 ) 

0178 


r\ r\ t , o r> 


GO TO oo 


o«> WPITT I 1 1 I * 2 0 ^ *3 ) 

r, 0 TO ^7 

oft wRITf < t:;.209ft) 

«>7 4R I TF t 1 1 ) * 2007 ) f . X t I ) 
09 «*RITF t ! U • 2U00 ) l')*NO 
MFtURN 
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?0Q0 

F 1 MAT ( 1 X/72H 

»«• 

L 

» 1 OR 

LESS./ 1 

7091 

FU iMAT ( 1X/22M 

»»« 

N 

• 0 OR 

l CSS*/) 

209* 

FORMAT I 1X/27H 

»•» 

IDENTICAL 

X VALUES./) 

209* 

FORMAT i 1X/13H 

*** 

X 

VALUfS 

OUT or SFQUENCT./) 


2097 

2009 


1 


F OP 4AT IftH 

FORMAT (6H 
36'< 


«♦ I 7 * 1 OX t ftMX ( I ) 

■ ♦ ! 7 * 1 OX * 1 HN - t I 7/ 


E12.3/ I 


ERROR DETLCTCD IN ROUTINL 


! N T P *o l / ) 


END 

FUNCTION BJ12.N) 


THIS ROUTINE COMPUTES The SPHERICAL BESSEL FUNCTIONS 
OF The 1ST KIND JNtZ) FOR Thf FIRST FIVE VALUES OF N 


NP 1 *N* 1 

GO TO < 10f20*30*l9*19» *NP1 
C 

C 0012 1 

10 IF ( 2 *NE *0 *0 I GO TO U 
BJ*1.0 

GO TO 60 

11 BJ*S I N I 2 ) / 2 
GO TO 60 

C 

19 ASSIGN 30 TO JN5WT 

GO TO 21 

C Jl < 2 ) 

20 ASSIGN 60 TO JNSWT 

21 IF I 2 «NF *0*0 *G0 TO 22 
BJ"0*0 

BJ1-0.0 
60 TO 23 

22 BJ»l$INtZ>- 2 *COS< 2 n/ 2**2 
BJ1-6J 

23 60 TO JNSWT»f60«3U) 

C J2I2I 

30 BJ«I l3.0-2*2 )*SIN(2 )-3,0»2*COS<2 ) >/Z«*3 

BJ2*BJ 

IFIN- 9 ) 60 * 40 • 40 

C 

C J3 < 2 I 
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